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ABSTRACT 

The  study  of  fractal  geometry  and  chaotic  dynamical 
systems  has  received  considerable  attention  in  the  past 
decade.   Motivated  by  the  interesting  computer  graphics 
produced  by  these  fields,  mathematicians  have  attempted  to 
formalize  the  theoretical  structure  of  the  results, 
physicists  have  attempted  to  apply  the  theory  to  real  world 
phenomena,  and  laymen  have  enjoyed  much  of  the  popular 
literature  and  television  programs  that  the  field  has 
fostered.   Unfortunately,  the  mathematics  associated  with 
these  subjects  has  made  them  inaccessible  to  most 
undergraduates,  even  if  they  have  a  strong  background  in 
mathematics.   This  thesis  presents  the  basic  ideas  of  fractal 
geometry  and  chaotic  dynamical  systems  in  a  setting  that  can 
be  understood  by  undergraduate  students  who  have  had  a  course 
in  advanced  calculus .   We  hope  it  will  allow  them  to  gain  an 
appreciation  of  the  fields  and  motivate  them  to  pursue 
further  study. 
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I.  INTRODUCTION 

The  subjects  of  fractals  and  chaos  have  attracted  considerable  attention  in 
the  last  decade.   This  interest  ranges  from  a  "cult  following"  of  laymen  who 
are  intrigued  by  the  intricate  computer  graphics  associated  with  the  fields,  to  a 
rigorous  mathematical  treatment  of  the  subjects  by  topologists  and  experts  in 
dynamical  systems,  and  to  applications  of  these  results  to  the  real  world  by 
engineers  and  physicists.  However,  these  subjects  have  been  almost  wholly 
inaccessible  to  undergraduates  because  of  the  level  of  mathematics  required  to 
study  them.  This  thesis  presents  the  subjects  of  fractals  and  chaos  in  a  setting 
that  can  be  understood  by  a  typical  undergraduate  student  with  a  solid 
background  in  mathematics  through  advanced  calculus. 

The  subjects  of  fractals  and  chaos  are  not  new.  The  German 
mathematician  Georg  Cantor  (1845-1918)  knew  about  fractals,  and  the  French 
mathematician  Jules  Henri  Poincare  (1854-1912)  knew  about  chaos  in 
dynamical  systems  in  the  late  nineteenth  century.   Additionally,  the  French 
mathematicians  Pierre  Fatou  (1878-1929)  and  Gaston  Julia  (1893-1978)  knew 
about  Julia  sets  in  the  1920s.  However,  it  was  not  until  the  1970s  that  high- 
speed computers  allowed  others  to  see  what  these  men  had  discovered  and  to 
recognize  the  true  potential  of  these  fields.  The  growth  of  these  fields  has 
resulted  in  significant  scientific  advances  in  the  past  decade.  While  the 
discovery  of  quantum  mechanics  and  relativity  had  a  profound  impact  on 
very  specialized  areas  of  science,  fractals  and  chaos  have  had  a  universal  effect 
on  the  whole  scientific  community.  Until  recently,  science  had  become  so 


specialized  that,  not  only  were  physicists  and  mathematicians  not 
communicating  with  each  other,  but  molecular  biologists  were  not 
communicating  with  population  biologists  either.   The  science  of  chaos  has 
served  to  bring  together  the  entire  scientific  community,  including  physicists, 
electrical  engineers,  mechanical  engineers,  biologists,  economists, 
astronomers,  meteorologists,  medical  scientists,  and  of  course, 
mathematicians. 

To  explain  briefly,  fractals  and  chaos  help  to  describe  the  universe. 
Natural  objects,  from  crystals,  plants  and  geological  formations  to  weather 
patterns  and  galaxies,  seem  to  have  a  fractal-like  structure  which  eludes 
description  by  traditional  geometric  means.  Additionally,  dynamical 
(changing)  systems  in  the  real  world,  from  turbulence  in  fluids  to  fluctuating 
economic  trends  and  unpredictable  weather  patterns  that  previously  defied 
description,  are  now  being  understood  through  the  mathematics  of  chaos. 
While  these  fields  are  still  in  their  infancy,  the  potential  they  have  already 
demonstrated  seems  very  promising  in  increasing  our  understanding  of  the 
physical  world. 

The  mathematics  required  to  understand  this  thesis  includes  basic  courses 
in  calculus  and  linear  algebra.  Additionally,  due  to  some  of  the  examples,  a 
familiarity  with  numerical  analysis  is  helpful,  but  not  required.  Most  results 
(specifically  in  the  fractals  portion)  are  presented  in  the  Euclidean  plane,  with 
a  brief  mention  of  more  abstract  spaces  where  applicable.  Mathematical 
proofs  that  can  be  understood  with  this  basic  background  are  presented  in 
their  entirety.   Results  requiring  more  advanced  mathematics  are  referenced. 
Where  applicable,  results  are  explained  and  interpreted  in  a  basic  setting.  The 


sections  on  Julia  sets  and  the  Mandelbrot  set  require  a  small  amount  of 
complex  variable  theory.  Finally,  only  discrete  dynamical  systems  are  treated 
(with  the  exception  of  the  section  on  the  Lorenz  equations),  so  a  knowledge  of 
differential  equations  is  not  required. 

None  of  the  material  in  this  thesis  constitutes  original  research.   Instead, 
it  is  a  synthesis  of  some  of  the  best  written  works  on  the  subjects  of  fractals 
and  chaos,  interpreted  and  presented  in  such  a  way  as  to  be  understood  by  a 
typical  undergraduate  with  a  basic  mathematical  background.  Hopefully,  it 
will  serve  to  stimulate  further  interest  and  deeper  study  in  these  exciting  new 
fields. 

Specific  references  are  listed  in  the  sections  to  which  they  apply,  and  all  of 
the  figures  and  examples  taken  directly  from  references  are  so  cited  as  they 
appear.  However,  we  note  that  the  two  primary  references  for  this  thesis  are 
Barnsley  (1988)  for  the  chapter  on  fractals,  and  Devaney  (1989)  for  the  chapter 
on  chaos.   Additionally,  the  references  used  for  general  information 
throughout  the  thesis  are  Ross  (1980)  for  advanced  calculus,  Anton  (1987)  for 
linear  algebra,  Churchill  (1990)  and  Boas  (1987)  for  complex  variables,  and 
Giordano/Weir  (1991)  and  Arnold  (1987)  for  differential  equations. 


II.        BASIC  CONCEPTS 

A.  INTRODUCTION 

There  are  several  recurring  mathematical  concepts  that  appear  in  the 
study  of  both  fractals  and  chaos.  An  understanding  of  these  basic  concepts 
simplifies  this  analysis.  Instead  of  presenting  them  as  they  arise,  we  next 
discuss  these  concepts  as  background  material.  They  are  presented  in  a  strict 
mathematical  setting  which  is  easily  adapted  to  their  further  treatment 
within  the  settings  of  fractals  and  chaotic  dynamical  systems.  This 
background  material  is  common  to  most  of  the  references  cited  (either  as 
presented  material  or  assumed  knowledge),  but  the  presentation  here  most 
closely  follows  those  of  Barnsley  (1988)  and  Devaney  (1989). 

B.  METRIC  SPACES 

We  first  introduce  the  concept  of  metric  spaces.  Although  the  notion  of  a 
metric  space  defined  below  may  seem  abstract  and  unfamiliar,  the  most 
common  examples  of  the  real  line  and  the  Euclidean  plane  are  encountered 
in  calculus.  Many  results  from  fractal  geometry  can  be  applied  to  general 
metric  spaces,  but  for  the  purpose  of  simplifying  them,  they  are  presented  in 
their  most  natural  or  familiar  settings  (such  as  the  real  line  or  the  Euclidean 
plane). 

A  space  is  simply  a  collection  of  objects  (a  set)  with  some  additional 
structure  imposed  on  it.  In  a  metric  space,  the  additional  structure  is  a  metric, 
or  distance  function,  which  relates  every  two  elements  x  and  y  in  the  set  to  a 


unique  nonnegative  real  number  d(x,  y),  called  the  distance  between  them. 
Any  distance  function  (or  metric)  must  satisfy  the  following  three  axioms  for 
all  x,  y,  z  in  the  space: 

1.  Symmetry:  d(x,  y)  =  d(y,  x); 

2.  Positive  Definiteness:  d(x,  y)  >  0,  and  d(x,  y)  =  0  <=>  x  =  y; 

3.  Triangle  Inequality:  d(x,  y)  <  d(x,  z)  +  d(z,  y). 

A  metric  space  is  denoted  by  (X,  d),  where  X  is  the  set  and  d  is  the 
particular  distance  function.  The  metric  space  (R,  dn)  is  the  real  line  R  with 
d(x,  y)  =  I  x  -  y  I .  The  space  (R2,  d^)  is  the  Euclidean  plane,  where  X  =  RXR 
and  d(x,  y)  =  [(xi  -  yi)2  +  (x2  -  V2)2]1/2.  Whenever  possible,  we  confine  our 
examples  to  the  real  line  or  the  Euclidean  plane  with  these  metrics. 

We  also  need  the  concept  of  closure,  a  thorough  treatment  of  which  can 
be  found  in  most  advanced  calculus  texts.  The  closure  of  a  set  is  the 
intersection  of  all  closed  sets  containing  that  set.  An  equivalent  definition, 
often  presented  as  a  theorem,  is  that  if  a  sequence  of  points  in  a  set  S 
converges  to  a  point  x  in  the  space,  then  x  is  in  the  closure  of  S.  Thus, 
every  set  is  a  subset  of  its  closure.  As  an  example,  consider  the  set  Q  of 
rational  numbers.  The  number  V2  is  not  in  Q,  but  the  sequence  of  rational 
numbers  1.4, 1.41, 1.414, 1.4142, . . .  converges  to  VI,  so  V2  is  in  the  closure  of 
Q  (in  fact,  every  real  number  is  in  the  closure  of  Q). 

A  concept  frequently  associated  with  closure  is  that  of  denseness.  A  set  S 
is  dense  in  a  metric  space  (X,  d)  if  for  every  point  x  e  X  and  all  e  >  0,  there 
exists  a  point  s  €  S  such  that  d(s,  x)  <  e.  For  example,  the  rational  numbers 
form  a  dense  subset  of  the  real  line.  Another  way  to  state  this  is  that  between 


any  two  points  on  the  real  line,  there  exists  a  rational  point.  If  a  set  S  is 
dense  in  the  space  X,  then  the  closure  of  S  is  X. 

Much  of  our  work  on  fractals  and  chaotic  dynamical  systems  involves  a 
topological  property  of  some  metric  spaces  called  completeness.  A  Cauchy 
sequence  is  a  sequence  in  which  the  terms  get  arbitrarily  close  to  each  other, 
i.e.,  for  any  e  >  0,  there  exists  a  number  N  such  that  d(x,v  xm)  <  e  for  all 
m,  n  >  N.  A  complete  metric  space  is  one  where  any  Cauchy  sequence  of 
points  in  the  space  converges  to  a  point  in  the  space.  An  intuitive  example  is 
to  again  consider  the  space  of  rational  numbers  and  the  sequence  1, 1.4, 1.41, 
1.414, 1.4142, . . . ,  which  converges  to  the  irrational  Vz  Clearly,  terms  in  this 
sequence  get  arbitrarily  close  to  each  other;  however,  the  limit  value  (V2)  is 
not  in  the  space  of  rational  numbers,  so  the  sequence  does  not  converge  to  a 
point  in  the  space.  This  is  to  say  that  there  are  "holes"  in  the  space  of  rational 
numbers,  hence  it  fails  to  be  a  complete  space.  The  real  line,  on  the  other 
hand,  is  complete,  since  all  Cauchy  sequences  of  real  numbers  do  converge  to 
a  real  number.  The  spaces  (R,  d)  and  (R2,  d£)/  are  examples  of  complete 
metric  spaces.  A  rigorous  development  of  the  concept  of  completeness  and 
formal  proofs  of  the  completeness  of  the  spaces  just  mentioned  can  be  found 
in  most  advanced  calculus  texts. 

C  ITERATED  FUNCTION  SYSTEMS 

Fractals  can  be  created,  and  chaotic  dynamical  systems  understood, 
through  an  analysis  of  iterated  function  systems.   A  simple  one-dimensional 
case  of  such  a  system  is  demonstrated  by  entering  a  starting  value  in  a 
calculator  and  repeatedly  activating  a  single  function  key  (for  example,  the  x2 


key).  The  result  of  executing  this  idea,  with  a  starting  value  of  2  and  the  x2 
key,  is  2, 22, 24, 28, 216, . . . .  In  this  example,  the  system  diverges  to  infinity 
and  demonstrates  a  regular  and  predictable  behavior.   Another  example  of  an 
iterated  function  system  on  a  calculator  is  to  use  the  Vx  key  and  start  with 
any  positive  number.  Regardless  of  the  initial  value,  this  system  will 
converge  to  1.  Unfortunately,  not  all  iterated  function  systems  are  quite  so 
simple.  The  study  of  chaotic  dynamical  systems  attempts  to  explain  irregular 
behavior  in  more  complicated  iterated  function  systems. 

A  function,  or  map,  is  a  rule  which  assigns  to  each  element  in  a  specified 
domain  a  unique  element  in  a  codomain.   The  usual  notation  for  this  concept 
is  f:  D-»R,  where  f  is  the  rule,  D  is  the  domain,  and  R  is  the  codomain.  A 
real-valued  function  of  a  real  variable  is  specified  by  f:  R-»R,  and  a  function 
mapping  the  closed  unit  interval  to  itself  is  written  f:  [0,  l]->[0, 1].  A  function 
f :  R-»R  of  the  form  f(x)  =  ax  where  a  is  a  constant,  is  said  to  be  linear;  if  it  is 
of  the  form  f (x)  =  ax  +  b  with  a  and  b  constants,  it  is  called  af  f  ine. 

A  function  is  one-to-one  if  each  element  in  the  range  is  the  image  of  a 
unique  element  in  the  domain;  that  is  f(x)  =  f(y)  implies  x  =  y.  A  function  is 
onto  if  every  element  in  the  codomain  is  the  image  of  at  least  one  element  in 
the  domain.  A  continuous  function  f:  R-»R  is  stated  as  f  e  C,  and  a 
continuously  differentiable  function  is  stated  as  f  e  C1.  A  function  that  is 
n-times  continuously  differentiable  is  expressed  by  f  €  O1,  and  its  nth 
derivative  is  written  as  f(n>(x).  If  a  function  is  one-to-one,  then  its  inverse  H 
exists  according  to  the  rule  fHx)  =  y  where  f(y)  =  x.  The  domain  of  f1  is  the 
range  of  f.   A  function  which  is  continuous,  one-to-one,  and  onto  for  which 
f-1  is  also  continuous  is  called  a  homeomorphism.  If  f  is  a 


homeomorphism  and  f  and  f-1  both  are  differentiable,  then  f  is  called  a 
diffeomorphism. 

As  previously  mentioned,  an  iterated  function  system  results  from 
repeatedly  applying  a  function  to  an  initial  point  x.  The  sequence  x,  f(x), 

f(f(x)),  f(f(f(x))), ...    is  also  written  as  x,  f(x),  f2(x),  f3(x), This  sequence  of 

iterates  is  called  the  forward  orbit  of  x  under  f,  and  is  denoted  by  0+(x). 
From  our  example  above,  the  forward  orbit  of  2  under  f(x)  =  x2  is  2, 4, 16, 

256, If  f  is  a  homeomorphism,  f1  exists,  and  the  backward  orbit  of  x 

under  f  is  the  sequence  x,  fHx),  f  2(x), . . . ,  denoted  by  O(x).  Since  f(x)  =  x2 
is  not  a  homeomorphism,  there  is  no  backward  orbit  under  f.  However,  if 

we  let   f(x)  =  2x,  then  the  backward  orbit  of  2  under  f  is  2, 1, 1/2, 1/4, A 

point  x  for  which  f(x)  =  x  is  called  a  fixed  point  of  f;  a  point  x  for  which 
f^Kx)  =  x  is  called  a  periodic  point  of  period  n.  The  smallest  positive  integer 
n  such  that  f  n(x)  =  x  is  called  the  prime  period  of  x.  A  point  x  is  said  to  be 
eventually  periodic  if  x  is  not  periodic,  but  there  exists  an  integer  m  such 
that  fnKx)  is  a  periodic  point.  For  a  function  f  where  the  derivative  is 
defined  on  the  entire  domain,  a  point  x  is  a  critical  point  if   f  (x)  ■  0,  for 
example  f(x)  =  x2  at  x  =  0.  It  is  a  degenerate  critical  point  if  it  is  a  critical  point 
and  f"(x)  =  0,  where  the  prime  notation  denotes  differentiation  of  the 
function. 

EXAMPLE.  The  following  simple  example  illustrates  a  situation  when  a  fixed 
point  is  guaranteed  to  exist  for  a  given  function.  Let  f:  [0,  l]-»[0, 1]  be  such 
that  f  €  C.  Then  f  has  at  least  one  fixed  point  on  [0, 1],  as  we  now  prove.  If 
f(0)  =  0  or  f(l)  =  1,  then  0  or  1  is  a  fixed  point  and  we  are  done.  Otherwise, 
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f(0)  >  0  and  f(l)  <  1,  and  define  g(x)  =  f(x)-x.  Since  g  is  the  difference  of  two 
continuous  functions,  g  is  continuous.  Also  g(0)  >  0  and  g(l)  <  0.  Hence,  by 
the  Intermediate  Value  Theorem  from  elementary  calculus,  there  exists  a 
point  a  e    [0, 1]  such  that  g(a)  =  0.  Therefore,  f(a)  =  a,  completing  the  proof 
(see  Figure  2.1).  This  proof  can  be  extended  to  any  continuous  function 
mapping  a  closed  interval  [a,  b]  to  itself. 


y  =  x 


0  a  1  X 

Figure  2.1  The  fixed  point  of  a  graph  f:  [0,  l]-»[0, 1]. 


We  will  frequently  use  a  geometric  or  graphical  analysis  to  investigate 
iterated  function  systems.  For  example,  the  orbit  of  two  points,  p  and  q, 
under  f(x)  =  x2  is  shown  graphically  in  Figure  2.2.  Clearly,  1  and  0  are  fixed 
points  of  f(x)  =  x2.  The  orbits  of  points  greater  than  1  diverge  to  infinity,  and 


the  orbits  of  points  in  the  interval  [0, 1)  converge  to  0.  This  example 
suggests  the  notion  of  fixed  points  (or  periodic  points)  as  being  "attracting"  or 
"repelling." 


f(x)  A 


y  =  x 


0  p      1q  X 

Figure  2.2  The  orbits  of  two  points  p  and  q  under  f(x)  =  x2. 


A  periodic  point  p  of  prime  period  n  (to  include  fixed  points  with  n  =  1) 
is  called  hyperbolic  if  I  (Wp)  1*1.  If  I  (Wp)  I  <  1,  then  p  is  called  an 
at  tract  or,  and  if  I  (WCp)  I  >1,  then  p  is  called  a  repellor.  As  in  the  above 
example  with  f(x)  =  x2,  the  orbit  of  points  near  an  attractor  tend  towards  that 
point,  while  the  orbits  of  points  near  a  repellor  tend  away  from  that  point.  A 
precise  definition  of  "near"  will  be  given  in  the  chapter  on  chaos.  With 
f(x)  =  x2,  the  derivative  is  f  (x)  =  2x,  so  the  derivative  evaluated  at  the  fixed 
points  0  and  1  gives  f'(0)  =  0  and  f  (1)  =  2.  Therefore,  0  is  an  attracting 
fixed  point  and  1  is  a  repelling  fixed  point.  Values  in  the  interval  [0, 1)  tend 


10 


towards  0  and  away  from  1,  while  values  in  the  interval  (1,  «>)  all  diverge. 
Hence,  the  point  0  "attracts"  iterates  and   1    "repels"  them. 

A  graphical  analysis  can  be  very  useful  when  analyzing  fixed,  periodic,  or 
eventually  periodic  points.   The  graphs  in  Figures  2.3  through  2.6  show 
functions  with  fixed  points,  periodic  points,  eventually  fixed  points,  and 
eventually  periodic  points.   While  geometric  constructions  cannot  replace 
rigorous  mathematical  proofs,  they  are  useful  to  demonstrate  and  provide 
insight  into  phenomena  that  occur  in  iterated  function  systems. 


y  =  x 


0  a  X 

Figure  2.3  A  fixed  point  a  of  f(x):  f(a)  =  a. 
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y  =  x 


Figure  2.4  Periodic  points  of  f(x):  f(pi)  =  p2  and  f(p2)  =  pi. 


y  =  x 


0      '     P 
Figure  2.5  An  eventually  fixed  point  p  of  f(x):  fCf^p))  =  f^p). 
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y  =  x 


f(p) 


Figure  2.6  An  eventually  periodic  point  p  of  f(x):  f2(p)  is  periodic. 

D.  CODE  SPACE 

A  useful  concept  for  identifying  different  points  on  fractals  (frequently 

called  addressing  points  on  fractals)   and  analyzing  chaotic  dynamical  systems 

is  that  of  code  space.  We  define  X2  as  the  set  of  all  infinite  sequences  of 

binary  digits  S1S2S3. . .  where  S{  e  {0, 1}.  We  next  define  a  distance  function 

for  all  s,  t  €  X2  by 

00 
d(s,t)  =  £lSi- til/2*. 
i=i 

Recalling  the  properties  of  a  metric,  it  is  easy  to  verify  that  this  mapping 

does  define  a  distance  function  on  X2  since   I  s\  - 1{  I  <  1.  .  Since  the  distance 

between  any  two  points  in  X2  is  dominated  by  the  convergent  geometric 

series 

i=i2 
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a  unique  and  finite  distance  exists  between  any  two  points  in  the  code  space 
(that  is,  d(s,  t)  is  well-defined). 

The  concept  of  code  space  can  be  extended  to  5jsj/  which  consists  of  the 
infinite  sequences  si,  S2,  S3, . . . ,  where  si  e  {0, 1, 2, 3, . . . ,  N-l).  Here  the 
distance  function  becomes 

ee 

d(s,t)=£lsi-til/(Ni). 

i  =  l 

It  is  an  easy  exercise  to  show  that  the  properties  of  a  metric  space  are  satisfied. 

We  will  frequently  perform  an  operation  on  code  space  known  as  the 
shift  map.  The  shift  map    a:  £2-^X2  is  given  by  a(siS2S3- . .)  =  (S2S3S4. . .).  The 
shift  map  simply  drops  the  first  term  in  the  sequence.  Since  g(0s2S3S4.  . .)  = 
a(ls2S3S4. .  .),  the  shift  map  fails  to  be  one-to-one. 

In  an  iterated  function  system  of  the  shift  map  a,  a1,  a2,  a3, . . .  on  Yji, 
the  only  fixed  points  are  000. . .  and  111. . . ;  the  eventually  fixed  points  are  of 
the  form  S1S2. .  .sn000. . .  and  siS2- .  .snlll. .  • ,  and  the  periodic  (or  eventually 
periodic)  points  consist  of  sequences  having  a  repeating  (or  eventually 
repeating)  block  of  0s  and  Is.  However,  there  are  infinitely  many  points 
which  are  neither  fixed  nor  periodic.  These  points  are  the  sequences  which 
have  no  repeating  nor  eventually  repeating  blocks  of  0s  and  Is. 

We  will  use  code  space  extensively  and  develop  it  in  greater  detail  when 
discussing  fractals  and  chaotic  dynamical  systems.  For  now,  this  introduction 
provides  a  brief  exposure  to  the  basic  ideas  of  this  important  concept. 
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E.    THE  CANTOR  SET 

The  Cantor  set,  which  has  traditionally  served  as  an  important 
pathological  example  in  analysis,  has  secured  a  more  distinct  role  in  the 
studies  of  chaos  and  fractals.  While  it  is  actually  one  of  the  simplest  examples 
of  a  one-dimensional  fractal,  we  first  present  it  in  its  traditional  construction 
to  ensure  familiarization  with  the  set  before  analyzing  it  in  its  new  guise. 

To  construct  the  classical  Cantor,  or  middle-thirds  set,  begin  with  the  unit 
interval  [0,1]  and  remove  its  middle  third  open  interval  (1/3,2/3).  Of  the 
two  remaining  line  segments,  [0, 1/3]  and  [2/3, 1],  remove  the  open  middle 
third  of  each  of  them.   Of  the  four  remaining  line  segments,  again  remove 
their  open  middle  thirds,  and  continue  this  iterative  procedure  ad  infinitum. 
The  remaining  set  of  points,  which  is  a  subset  of  the  unit  interval,  is  called 
the  Cantor  middle-thirds  (or  ternary)  set.  It  is  sometimes  referred  to  as 
"Cantor  dust"  because  of  the  scattered  configuration  of  the  remaining  points. 
A  picture  of  the  construction  of  the  classical  Cantor  set  is  shown  in  Figure  2.7. 

The  Cantor  set  clearly  contains  an  infinite  number  of  points  because,  at 
the  very  least,  the  endpoints  of  the  line  segments  left  after  each  iteration 
remain  in  the  set.   Nevertheless,  the  total  amount  of  length  removed  from 
the  unit  interval  is  equal  to  one.  To  see  this,  consider  the  length  removed 
during  each  iteration.  In  the  first  step,  a  segment  of  length  1/3  is  removed; 
in  the  second  step,  two  segments  each  of  length  1/9  are  removed;  and  in  the 
third  step,  four  segments  each  of  length  1  /27  are  removed.  More  generally, 
in  the  kth  step,  2k_1  segments  each  of  length  l/3k  are  removed.  Summing 
the  total  lengths  removed  yields 
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1/3  2/3 


Figure  2.7  Constructing  the  classical  Cantor  set  by  removing  the  middle 
third  through  infinite  iteration.  The  Cantor  set  consists  of  the  points  that 
remain. 


L  =  X  2n"1(l/3)n=  (l/3)£  (2/3)n=  (1/3)[1/(1  -  2/3)]  =  1. 

n=l  n=0 

So  the  total  length  removed  is  the  entire  length  of  the  unit  interval. 

Another  way  of  approaching  the  Cantor  set  is  through  ternary  expansion. 
Consider  expanding  each  point  in  the  unit  interval  in  its  base  three  form 
O.X1X2X3. . . ,  where   xi  e  {0, 1, 2}.  The  value  of  each  point  is  then  xi(l/3)  + 
X2(l/32)  +  X3O/33)  +. . . .  There  is  a  minor  technicality  with  this  approach, 
because  many  points  in  the  unit  interval  have  dual  ternary  representations. 

For  example,  the  point  1/3  can  be  expressed  both  as  .1000. . .  and  as  .0222 

In  fact,  the  endpoints  of  every  interval  in  our  previous  construction  of  the 
Cantor  set,  except  for  0  and  1,  have  such  a  dual  representation.  The  way  we 
remedy  this  is  to  always  use  the  representation  containing  the  repeating  2. 
Now,  all  of  the  points  between  1/3  and  2/3  have  xi  =  1,  and  these  are  the 
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points  removed  in  the  first  step  of  our  former  construction.   Similarly,  the 
points  removed  in  the  second  step  have  the  form  xi  =  0  or  2,  and  X2  =  1. 
Continuing  with  this  analysis,  we  see  that  the  points  remaining  in  the  Cantor 
set  are  precisely  those  that  have  no  Is  appearing  in  their  ternary  expansion. 

Thus  the  Cantor  set  consists  of  all  points  in  the  unit  interval  having 
only  Os  or  2s  in  their  ternary  expansion.  Moreover,  every  point  which  has 
only  Os  or  2s  in  its  ternary  expansion  is  a  member  of  the  Cantor  set.  This 
definition  proves  to  be  very  useful  when  we  analyze  chaotic  dynamical 
systems.  Note  that  the  Cantor  set  contains  no  interval  subsets  because 
between  any  two  points,  there  is  a  point  with  a  1  in  its  ternary  expansion.  To 
see  this,  assume  that  the  Cantor  set  does  contain  an  interval  subset.  We  can 
let  the  endpoints  of  this  interval  be  xiX2- .  -xn0xn+2. . .  and  xiX2. .  .xn2xn+2. . . , 
where  the  xi  agree  for  the  first  n  digits,  are  different  for  the  n+lst  digit 
(which  means  one  must  be  a  0  and  the  other  must  be  a  2),  and  after  which 
the  digits  can  be  arbitrary.  Since  every  point  between  these  two  must  be  in  the 
Cantor  set,  points  of  the  form  xiX2- .  .xnl. . .  must  be  in  the  set,  contradicting 
the  fact  that  the  Cantor  set  contains  no  points  with  Is  in  their  ternary 
expansion.  This  proves  that  the  Cantor  set  contains  no  intervals,  thus  its 
only  "connected"  subsets  are  single  points.  (We  discuss  the  concept  of 
connectedness  in  detail  in  the  chapter  on  Julia  sets).   This  illustrates  another 
property  of  the  Cantor  set:  it  is  totally  disconnected. 

Now  a  great  paradox  results.  We  have  removed  the  entire  length  of  the 
interval.  The  only  points  remaining  are  the  points  of  the  form  x  =  .X1X2X3  . . . 
where  xi  e  {0,  2}.  Now  we  form  the  function  f(x)  =  y  where  yi  =  xj/2.  The 
set  of  y  values  is  the  set  of  all  strings  of  Os  and  Is.  The  cardinality  of  this  set 


17 


is  the  same  as  the  cardinality  of  [0, 1]   since  there  is  a  one-to-one 
correspondence  with  the  binary  expansion  of  the  points  of  the  unit  interval! 
Thus,  the  cardinality  of  a  set  does  not  tell  the  whole  story  of  "how  many" 
points  are  actually  there  from  the  point  of  view  of  length  (in  the  one- 
dimensional  case). 

There  are  other  versions  of  Cantor  sets  resulting  from  similar 
constructions.   For  instance,  remove  instead  the  middle  fourth,  fifth,  or  some 
other  fixed  (1/N)  length  at  each  step.  We  can  also  remove  from  the  unit 
interval  an  open  middle  segment  of  length  a/3k  at  the  kth  stage  where 
0  <  a  <  1  is  fixed.  Then  we  can  show  the  total  length  removed  is  a.  These 
are  frequently  called  "fat"  Cantor  sets,  because  the  total  length  removed  from 
the  unit  interval  is  now  less  than  one.  These  sets  can  also  be  described  in 
terms  of  a  base  N  expansion,  from  which  it  can  be  shown  that  a  fat  Cantor  set 
contains  no  intervals. 
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III.  FRACTALS 

A.    THE  SETTING  FOR  FRACTALS 

Fractal  geometry  involves  the  study  of  certain  subsets  of  metric  spaces.  It 
can  be  viewed  as  an  extension  of  Euclidean  geometry  and  is  frequently  used  to 
describe  objects  that  occur  in  nature,  such  as  crystals,  plants,  clouds,  and 
geological  formations  (see,  for  example,  Cherbit  (1991)).  Fractals  have  been 
applied  to  computer  graphics  to  store  information  efficiently,  examples  of 
which  can  be  found  in  Barnsley  (1988).  Moreover,  the  use  of  fractals  to  study 
real-world  phenomena  has  provided  a  new  way  of  analyzing  the  world. 
Finally,  as  we  shall  see  in  Chapter  IV,  a  primary  use  of  fractals  is  to  classify 
and  analyze  chaotic  dynamical  systems. 

The  primary  references  for  this  chapter  are  Barnsley  (1988)  and  Falconer 
(1990).  While  the  presentation  most  closely  follows  that  of  Barnsley,  a  more 
rigorous  mathematical  development  of  most  of  the  results  can  be  found  in 
Falconer.  Additionally,  the  article  by  Harrison  (1989)  covers  much  of  the 
same  material,  while  many  of  the  examples  presented  here  are  from  Cherbit 
(1991). 

Much  of  the  current  literature  differs  in  the  precise  definition  of  a  fractal, 
so  our  approach  is  to  develop  the  setting  for  the  space  in  which  fractals  exist. 
Then  we  provide  many  examples  of  fractals  in  that  setting.  This  approach 
provides  a  good  initial  understanding  of  fractals  without  the  expense 
involved  in  achieving  a  thorough  understanding  of  their  elusive  definition. 
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A  useful,  although  incomplete,  definition  of  a  fractal  is  that  "a  fractal  is  a 
fixed  point  of  a  certain  kind  of  transformation  on  the  space  (H(X),  h(d))." 
Thus,  we  must  first  define  the  set  H(X)  and  the  distance  function  h(d). 
(H(X),  h(d))  is  a  metric  space  obtained  from  a  complete  metric  space  (X,  d), 
where  the  points  in  H(X)  are  closed,  bounded,  nonempty  subsets  of  X,  and 
where  h(d)  is  a  distance  function  based  on  the  metric  d,  (which  we  define 
shortly).  The  points  in  H(X)  are  called  the  "compact"  subsets  of  X.  To  avoid 
the  necessity  of  employing  concepts  from  advanced  calculus,  we  consider  only 
Euclidean  spaces  where  the  compact  sets  are  the  closed  and  bounded  subsets. 

Consider  the  Euclidean  plane  (R2,  dfi).  We  denote  by  H(R2)  the 
collection  of  all  closed  and  bounded  subsets  of  R2,  excluding  the  empty  set. 
Hence,  the  closed  unit  square  [0, 1]X[0, 1]  belongs  to  K(R2),  as  does  the  origin; 
however,  the  interval  (0, 1)X{0)  is  not  an  element  in  H(R2)  since  it  is  not  a 
closed  subset  of  R2.  Clearly,  the  union  of  any  two  elements  in  H(R2)  again 
belongs  to  H(R2).  The  intersection  of  two  elements  in  §i(R2)  is  not 
necessarily  an  element  of  Ji(R2)  as  the  intersection  may  be  empty. 

In  order  to  create  a  metric  space  out  of  the  set  K(R2),  a  distance  function 
(metric)  is  required  that  relates  any  two  elements  of  H(R2)  to  a  nonnegative 
real  number.  To  this  end  the  concept  of  dilation  is  helpful.  Given  any  closed 
subset  A  of  R2  and  e  >  0,  the  e-dilation  of  A  is  defined  to  be  the  set  of  all 
points  x  in  R2  such  that  the  smallest  Euclidean  distance  between  x  and  any 
point  in  the  set  A  is  less  than  or  equal  to  e;  i.e.,  the  e-dilation  of  A  is  the  set 
{x:  dE(x,  A)  <  e}.  For  example,  the  e-dilation  of  the  origin  is  simply  the  closed 
disc  with  radius  e.  The  e-dilation  of  the  unit  square  [0, 1]X[0, 1]  with    e  =  1/2 
is  shown  in  Figure  3.1.  The  points  of  the  square  also  belong  to  its  e-dilation. 
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Figure  3.1  The  e-dilation  of  the  unit  square  with  e  =  1/2. 

We  now  define  the  distance  from  an  element  A  in  H(R2)  to  an  element 
B  in  H(R2)  as  the  number  dh(A,  B)  representing  the  smallest  e  such  that 
every  point  in  B  is  covered  by  the  e-dilation  of  the  element  A.  In  other 
words,  dh(A,  B)  =  min(e:  y  e  {the  e-dilation  of  A}  V  y  e  B).  Clearly,  dh 
exhibits  the  second  and  third  properties  of  a  metric  from  the  original 
definition.  However,  as  shown  in  Figure  3.2,  it  is  not  symmetric  since 
dh(A,  B)  *  dh(B,  A)  in  general.  To  remedy  this  difficulty,  we  define  h(A,  B)  = 
max{dh(A,  B),  dh(B,  A)}.  The  number    h(A,  B)  does  satisfy  the  properties  of  a 
metric  for  points  in  the  set  H(R2),  and  is  referred  to  as  the  Hausdorf f  distance. 
Now  the  set  H(R2)  together  with  the  metric  h(A,  B)  is  a  metric  space.  It  is 
denoted  by  (H(R2),  h(d)),  or  for  purposes  of  simplicity  just  ft,   and  we  refer  to 
it  as  the  Hausdorf  f -metric  space.  We  remark  that  if  the  metric  space  (X,  d)  is 
complete,  then  the  associated  Hausdorff-metric  space  (K(X),  h(d»  is  also 
complete.  In  particular,  H  is  complete  since  X  =  R2  is  complete. 
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dh(B,  A) 


Figure  3.2  The  distances  dh(A,B)  and  dh(B,A)  are  not  equal  in  general. 

B.    CONTRACTION  MAPPINGS 

In  order  to  generate  fractals  geometrically  we  need  to  be  familiar  with  the 
concept  of  contraction  mappings.    The  results  of  this  section  are  presented  in 
the  setting  of  a  general  metric  space.  However,  we  concentrate  our  efforts  on 
three  main  spaces:  the  real  line  with  the  normal  distance  function,  the 
Euclidean  plane,  and  the  Hausdorff-metric  space  (H(R2),  h(dE». 

A  mapping  f:  X— »X   is  said  to  be  a  contraction  mapping  if  there  exists  a 
constant   0  <  s  <  1  such  that  d(f(x),  f(y»  <  s(d(x,  y))  for  all  x,  y  e  X.  The 
number  s  is  called  the  contractivity  factor  for  f.  We  prove  shortly  that  a 
contraction  mapping  is  always  continuous.  Under  a  contraction,  any  two 
points  in  the  space  that  begin  a  distance  D  apart  will  be  moved  to  within  a 
distance  sD  of  each  other.  A  key  result  concerning  contraction  mappings, 
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and  one  which  is  critical  when  constructing  fractals  as  subsets  of  the 
Hausdorff-metric  space  Ji,  is  discussed  next. 

THE  CONTRACTION  MAPPING  THEOREM. 

If   f:  X— >X  is  a  contraction  defined  on  a  complete  metric  space  (X,  d),  then  f 
has  a  unique  fixed  point  xf  €  X.  Furthermore,  for  any  point  x  €  X,  the 
sequence  f°(x)  converges  to  Xf. 

To  prove  this  we  first  need  the  following  results. 

THEOREM.  If  f  is  a  contraction,  then  it  is  continuous. 
PROOF.  If  we  let  s  be  the  contractivity  factor  of  f,  and  if  e  >  0  is  given,  then 
choosing  5  =  e/s  yields  d(x,  y)  <  8  =>  d(f(x),  f(y))  £  s(d(x,  y))  <  s(e/s)  =  e.  This 
shows  that  f  is  continuous. 

LEMMA.  If   f   is  a  contraction  with  contractivity  factor  s,  then  for  a  fixed  x 
and   m  <  n,   dWx),  f^Cx))  <,  smdCx,  f*-m(x)). 

PROOF.  The  proof  follows  immediately  from  the  contractivity  factor  s,  and 
the  principle  of  mathematical  induction. 

PROOF  OF  THE  CONTRACTION  MAPPING  THEOREM. 

Let  xo  be  an  arbitrary  point  in  the  complete  metric  space  (X,  d),  and  let 
f:  X-»X  be  a  contraction  mapping  such  that  xi  =  f(xo),  X2  =  f(xi)  =  f^Cxo),  and  in 
general,  xn  =  f(xn.i)  =  f^Uo).  For  positive  integers  m  and  n  such  that  m  <  n, 
we  have 
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d(xnv  xn)  =  d(f^(xo),  f^(xo))  =  dCf^Xxo),  W^xo))) 

<  smdCxo,  ^-"Kxo))  =  S^xo,  xn_m) 

<  sm[d(xo,  xi)  +  d(xi,  X2)  +. . .+  d(xn.m-i,  xn.m)] 

(by  the  triangle  inequality  of  d) 

<  smd(xo,  xi)[l  +  s  +  s2+. . .+  s^™-1] 

(by  applying  the  lemma  to  each  term) 

<  [sm/(l  -  s)]d(xo,  xi)     (by  the  geometric  series). 

Since  s  <  1,  s™-*  0,  hence  d(xIIV  xn)-»0.  Since  X  is  a  complete  metric  space, 
xn-»xf  for  some  xf  e  X,  so  the  sequence  (f^Kxo)}  converges  to  the  point  xf  for 
any  xo  in  the  space.  To  show  that  xf  is  a  fixed  point    since  f  is  a  continuous 
function,  f(xf)  =  f(lim  xn)  =  lim  f(xn)  =  lim  xn+i  =  xf.    Therefore,  xf   is  a  fixed 
point  of  f . 

To  show  that  xf  is  the  only  fixed  point  of  f,  assume  that  y  is  also  a  fixed 
point.  Then  d(xf,  y)  -  d(f(xf),  f(y))  <  sd(xf,  y).  Since  s  <  1  and  d  £  0,  this 
implies  that  d(xf,  y)  =  0,  so  xf  =  y,  completing  the  proof. 

The  importance  of  the  contraction  mapping  theorem  will  become 
apparent  when  fractals  are  created  through  iterated  function  systems.  When 
we  find  a  point  in  the  space  that  is  a  fixed  point  under  a  contraction  mapping, 
we  then  know  it  is  unique  (hence  we  will  refer  to  "the"  fixed  point  of  the 
mapping).   Furthermore,  we  know  that  any  initial  point  converges  to  the 
fixed  point  under  iteration  of  f. 
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C    AFFINE  TRANSFORMATIONS  OF  THE  PLANE 

Affine  mappings  were  defined  earlier  as  functions  of  the  form 
f(x)  =  ax  +  b.  In  R2,  we  change  our  notation  slightly  and  write  an  affine 
transformation  as  w(x)  =  Ax  +  t,  where  x  is  a  point  in  R2  (a  column  two- 
vector),  A  is  a  2X2  matrix,  and  t  is  a  (fixed)  vector  in  R2.  We  frequently 
use  the  form 

w(> 


--GH:  *;]#" 


The  mapping  w(x)  =  Ax  is  called  a  linear  transformation,  and  it  takes  any 
parallelogram  with  one  vertex  at  the  origin  into  another  parallelogram  with  a 
vertex  at  the  origin,  provided  that  the  determinant  of  A  is  not  zero  (that  is, 
A  is  nonsingular).  Notice  that  the  origin  remains  fixed  under  any  linear 
transformation.  A  result  from  linear  algebra  is  that  the  matrix  A  can  be 
written  in  the  form 


[:!]■ 


rjcosdj     -r2sin82 
rjsinGj      r2pos92 


(1) 


where  the  point  (a,  c)  has  polar  coordinates  (ri,  0i)  and  the  point  (b,  d)  has 
polar  coordinates  fo,  82  +  n/2).  The  transformation  w(x)  =  Ax  "deforms" 
the  space  R2  relative  to  the  origin.  By  deform,  we  refer  to  the  result  from 
linear  algebra  that  any  linear  transformation  can  be  expressed  as  a  succession 
of  shears,  compressions,  expansions,  and  rotations  of  the  space.  The  result  of 
adding  the  vector  t  to  a  point  in  R2  shifts  that  point  by  the  magnitude  and 
direction  of  t,  so  an  affine  transformation  of  the  plane  can  be  thought  of  as  a 
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deformation  (shearing,  shrinking /stretching  and  rotating)  followed  by  a 
translation.   Figure  3.3  shows  the  result  of  the  affine  transformation 


■EH ° -©>[■] 


on  the  triangle  T  with  vertices  A  =  (0, 1),  B  =  (1, 1),  and  C  =  (1,  0),  where  A', 
B',  and  C  are  the  images  of  A,  B,  and  C  under  w,  respectively. 


1/2 


0     x 


Figure  3.3  An  affine  transformation  of  a  triangle. 


We  now  consider  several  special  affine  transformations.   From  the 
representation  in  polar  coordinates  Eq.  (1)  above,  if  ri  =  T2  =  r,  and  Q\  =  6  2=  ®* 
we  say  the  transformation  is  a  similitude,  and  write  it  in  the  form 


or, 


-M- 


rcos6  -rsinG 
rsinG  rcos6j 


][:;]*» 


{:;]■ 


rcos8  rsin8 
rsin6  -rcosG 


fcPM 
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A  similitude  scales  any  image  in  R2  by  the  same  factor  r  in  each  component 
direction,  rotates  it  by  the  angle  0,  and  then  translates  it  by  the  magnitude  and 
direction  of  the  vector  t.   Recall  from  linear  algebra  that  the  transformation 


cos0  -sinG 
sinG  cosG 


is  a  rotation,  and  the  transformation 


[-] 


is  a  reflection  across  the  x-axis  in  R2.   Both  of  these  transformations  behave  as 
one  would  expect,  based  on  their  names.  Figure  3.4  shows  the  transformation 


w 


[:;]■ 


cos(tc/2)      -sin(7c/2) 
sin(rc/2)      cos(tc/2) 


M 


(2) 


of  the  unit  square  which  is  a  rotation  with  0  =  jc/2.  Figure  3.5  shows  the 
transformation 


w 


[:;]-[;  ■"][:;] 


(3) 


of  the  unit  square  which  is  a  reflection  about  the  x  axis.  Being  familiar  with 
both  reflections  and  rotations  is  very  useful  when  creating  fractals. 
The  action  of  w(x)  =  Ax  +  t  can  be  determined  by  evaluating  its  effect  on  any 
point  (xi,  X2)  in  the  plane.  However,  it  is  sometimes  more  useful  to 
determine  which  affine  transformation  has  caused  a  particular  (observed) 
change  to  a  given  set  of  points.  The  beauty  of  affine  transformations  is  that 
they  are  uniquely  determined  by  their  action  on  any  three  non-collinear 
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0    x 


1 


Figure  3.4  A  rotation  of  the  unit  square  by  the  transformation  in  Eq.  (2). 
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Figure  3.5  A  reflection  of  the  unit  square  by  the  transformation  in  Eq.  (3). 

points  in  the  plane.  To  see  this,  start  with  any  three  points  (xi,  X2>,  (yi,  y2)/ 
and  (zi,  Z2>   and  note  their  movement  under  an  (unknown)  affine 
transformation  to  the  new  points  (xi*,  X2'),  (yi',  yi'),  and  (zi*,  Z2')«  The 
coefficients  a,  b,  c,  d,  e,  and  f  of  the  transformation  can  be  determined  from 
the  system  of  linear  equations 

xia  +  X2b  +  e  =  xi', 

yia  +  y2b  +  e  =  yi', 

zia  +  Z2b  +  e  =  zi', 
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xjc  +  X2d  +  f  =  X21, 
yic  +  y2d  +  f  =  yi, 

7.\Q  +  Z2d  +  f  =  Z2'. 

Because  the  points  are  non-collinear,  the  matrix  of  the  coefficients  is 
nonsingular  and  the  system  has  a  unique  solution.   An  example  follows. 

EXAMPLE.  Using  the  example  in  Figure  3.3,  the  transformation  on  the 
triangle  transforming  the  points  (1, 0),  (1, 1),  and  (0, 1)  to  (-2, 1/2),  (0, 1/2), 
and  (0, 1)  respectively  yields  the  following  system  of  equations: 

(0)a  +  (l)b  +  e  =  -2, 

(l)a  +  (l)b  +  e  =  0, 

(l)a  +  (0)b  +  e  =  0, 

(0)c  +  (l)d  +  f  =  l/2, 

(l)c  +  (l)d  +  f  =  l/2, 

(l)c  +  (0)d  +  f  =  l. 
Solving  for  the  coefficients  a  through  f  yields  a  =  2,  b  =  0,  c  =  0,  d  =  -1/2, 
e  =  -2,  and  f  =  1,  hence,  the  affine  transformation  is  given  by 


{:;]■[»  -EH?} 


D.    CONTRACTION  MAPPINGS  OF  THE  SPACE  K(X) 

Given  n  continuous  maps  of  a  metric  space  wi,  W2, . .  • ,  wn:  X-»X,  we 
construct  a  map  W  of  the  associated  space  K   by  considering 

W(x)  =  wi(x)  u  W2(x)  u. .  .u  wn(x), 
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where  x  is  any  point  in  the  space  H,  i.e.,  x  is  a  closed,  bounded,  nonempty 
subset  of  the  original  metric  space.  For  x  e  H,  we  define  wi(x)  =  {wi(y):  y  e  x 
c  X}.  We  state  here  that  the  mapping  W  maps  H   to  itself.  This  fact  is  based 
on  two  important  results  from  advanced  calculus,  which  we  state  but  do  not 
prove.   The  first  result  is  that  continuous  images  of  compact  sets  are  compact, 
so  that  each  wj  maps  H  to  itself.  The  second  is  that  the  finite  union  of 
compact  sets  is  compact.  Thus,  W  maps  H  to  itself,  as  claimed.  For  a  simple 
illustration  of  this  concept,  let  n  =  2  and  X  =  R,  and  consider  the  maps  on  R 
defined  by    wi(x)  =  x  and  W2(x)  =  2x.  If  we  start  with  the  singleton  compact 
set  x  =  {1},  we  have  W(l)  =  wi({l})  u  w2({l})  -  {1}  u  {2}  =  {1, 2}. 

We  next  iterate  the  map  W.  The  set  of  points  Wn(x)  grows  quite  rapidly 
as  n  increases.  To  see  this,  consider  the  same  example  used  above  with  n  =  2. 
Then  W2(x)  =  W(W(x))  =  W(wi(x)  u  W2(x))  =  wi(wi(x)  u  W2(x))  u  W2<wi(x)  u 
W2(x))  =  wi(wi(x))  u  wi(w2(x))  u  W2(wi(x))  u  W2(w2(x)).  Using  wi  and  W2 
as  above,  and  again  starting  with  x  ■  {1},  we  have  W2({1})  =  W({1,  2})  = 
{1, 2}  u  {2, 4}  =  {1, 2, 4).  Similarly,  W3({1})  =  {1, 2, 4}  u  {2, 4, 8}  =  {1, 2, 4, 8). 
Note  that  wi(w2(x))  =  W2(wi(x))  because  wi  and  W2  are  linear 
transformations  on  R,  which  commute.   For  notational  convenience,  we 
shall  not  distinguish  between  the  notation  x  ■  1  and  x  =  {1}  for  the  space 
K(R)  includes  sets  of  singleton  points. 

We  use  the  abbreviated  notation  Wi(wj(x))  =  wij(x)  etc.  Then, 
W1  =  wi  u  W2;  W2  =  wii  u  W12  u  W21  u  W22/  and  W3  =  win  u  wn2  u  W121  u 
W122  vj  W2ii  u  W212  u  W221  u  W222-  These  subscripts  may  look  familiar.  We 
are  building,  in  this  way,  an  iterated  function  system  in  one-to-one 
correspondence  with  the  code  space  Yji-  To  make  this  clear,  we  let  each  1  in 
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our  subscripts  correspond  to  the  0  in  code  space,  and  let  each  2  in  our 
subscripts  correspond  to  the  1  in  code  space.  Then,  the  iterate  Wijk...  is  in 
direct  correspondence  with  the  point  i-1,  j-1,  k-1, ...  in  £2-  We  could  have 
emphasized  this  more  had  we  named  our  contractions  wo  and  Wj. 
However,  we  save  the  designation  wo  for  a  specific  use  later  on.  Similarly,  if 
W  is  the  union  of  n  mappings,  then  the  subscripts  are  analogous  to  the 
points  in  £n  after  infinite  iteration.  This  concept  proves  quite  useful  when 
we  identify  points  on  fractals  through  an  addressing  scheme. 

Notice  that  the  infinite  iteration  W°°(l)   in  our  simple  example  produces 
the  set  {1,  2,  4,  8, 16, . . .}.  While  the  dynamics  of  this  particular  iterated 
function  system  are  not  very  exciting,  when  each  wi  is  a  contraction  with 
associated  contractivity  factor  Si  the  results  do  become  much  more 
interesting,  as  we  discuss  next. 

To  discuss  contraction  mappings  on  the  space  H,  we  need  the  following 
two  theorems. 

THEOREM.  If  a  mapping  on  a  metric  space  w:  X-»X  is  a  contraction  with 
contractivity  factor  s,  then  w:  Jt(X)-»Jt(X)  is  also  a  contraction  with  contrac- 
tivity factor  s. 
PROOF.  Let  B,  C  e  Ji.  Then 

dh(w(B),  w(C»  =  min{e:  y  e  {the  e-dilation  of  w(B)}  V  y  €  w(C)} 
<  min{se:  y  e  {the  e-dilation  of  B}  V  y  e  C} 
=  s  min{e:  y  e  {the  e-dilation  of  B}  V  y  €  C) 
=  s  dh(B,  C). 
Similarly,  dh(w(C),  w(B»  <  s  dh(C,  B). 
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Therefore,  h(w(B),  w(C))  <  s  h(B,  C),  completing  the  proof. 

THEOREM.  Let  w\,  w2, . . . ,  wn  be  contractions  on  the  space  H  with 

contractivity  factors  si,  S2, . . . ,  sn.  Define  W:  H-tft   by 

W(B)  =  wi(B)  u  W2(B)  u  . . .  u  wn(B) 

for  each  B  e  H.  Then  W  is  a  contraction  mapping  on  Ji  with  contractivity 

factor  s  =  max{si:  i  =  1, 2, . . . , }. 

PROOF.  We  prove  this  for  the  case  n  =  2.  The  rest  follows  by  induction. 

Let  B,Ce  H.  Then 

H(W(B),  W(C))  =  h(wi(B)  u  w2(B),  wi(C)  u  w2(C)) 

£  max{h(wi(B),  wi(C)),  h(w2(B)/  w2(C))} 
<,  max{sih(B,  C),  s2h(B,  C)}  =  sh(B,  C). 

A  simple  exercise  shows  why  the  first  inequality  holds.  This  completes  the 

proof. 

Figure  3.6  shows  two  contractions  on  the  unit  square  in  R27  wi  contracts  by 
1/2  in  the  xi  direction,  and  w2  contracts  by  1/4  in  the  x2  direction.  Their 
union  has  a  contractivity  factor  of  1/2  by  the  previous  theorem. 

Since  W  is  a  contraction  on  the  complete  metric  space  K(X)  we  know 
that  it  has  a  unique  fixed  point  in  Jt(X),  which  is  a  closed,  bounded,  non- 
empty subset  of  the  space  X.  By  the  Contraction  Mapping  Theorem,  all  points 
in  K  tend  to  this  fixed  point  under  infinite  iteration  of  W.  This  unique 
fixed  point  is  called  the  attractor  of  W,  and  frequently  exhibits  very  interesting 
properties. 
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Figure  3.6  The  union  of  two  contractions  of  the  unit  square. 

E    CREATING  FRACTALS  THROUGH  ITERATED  FUNCTION  SYSTEMS 
1.      The  Cantor  Set 

We  are  now  prepared  to  construct  fractals  using  affine  transfor- 
mations and  the  Contraction  Mapping  Theorem.  We  start  with  the  simplest 
setting,  namely  the  unit  interval,  and  construct  a  fractal  with  which  we  are 
already  familiar. 

Consider  the  union  of  two  affine  transformations  of  the  real  line  R, 
given  by  W(x)  =  wi(x)  u  W2(x),  where  wi(x)  =  x/3,  and  W2(x)  =  x/3  +  2/3.  If 
we  first  restrict  our  attention  to  the  unit  interval  [0, 1],  then  under  the  first 
iteration  of  wi  the  unit  interval  is  shrunk  by  1/3  towards  the  origin  such 
that  [0,  l]-»[0, 1/3].  Under  the  first  iteration  of  W2  the  unit  interval  is  shrunk 
by  1/3  and  then  translated  to  the  right  by  2/3  so  that  [0,  l]-»[2/3, 1].  Note 
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that  the  effect  of  W  on  [0, 1]  is  the  removal  of  the  open  middle  third  of  the 
unit  interval,  the  first  step  in  constructing  the  Cantor  "middle  thirds"  set. 

Under  the  second  iteration  W2,  the  mappings  wi  and  w2  now  act  on 
the  subsets  [0, 1/3]  and  [2/3, 1]  of  the  unit  interval,  so  that  W2[0, 1]  =  wn[0, 1]  u 
wi2[0, 1]  u  W2i[0, 1]  u  W22[0, 1].  This  can  be  written  as  wi[0, 1/3]  u  wi[2/3, 1]  u 
w2[0, 1/3]  u  w2[2/3, 1].  Since  wi[0, 1/3]  =  [0, 1/9],  wi[2/3, 1]  =  [2/9, 1/3], 
w2[0, 1/3]  =  [2/3, 7/9],  and  w2[2/3, 1]  =  [8/9, 1],  we  obtain  VV^O, 1]  = 
[0, 1/9]  u  [2/9, 1/3]  u  [2/3, 7/9]  u  [8/9, 1].  This  procedure  has  the  same  effect  as 
removing  the  open  middle  third  of  each  interval  [0, 1/3]  and  [2/3, 1]. 

Using  a  similar  analysis,  W3[0, 1]  again  removes  the  open  middle 
third  of  each  of  the  four  intervals  produced  from  W2,  and  continuing  in  this 
manner  results  in  the  construction  of  the  Cantor  middle-thirds  set  (see  Figure 
3.7). 
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Figure  3.7  Creating  the  Cantor  set  through  an  iterated  function  system  (IFS). 

The  above  construction  is  perhaps  the  simplest  example  of  a  fractal. 
Further  examination  reveals  that  this  set  has  some  interesting  properties. 
First,  consider  the  set  of  points  that  are  fixed  under  the  map  W.  Note  that  the 
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fixed  point  of  wi  is  0,  since  (l/3)(0)  =  0,  and  that  the  fixed  point  of  W2  is  1, 
since  (1/3)(1)  +  2/3  =  1.  Since  each  of  these  mappings  is  a  contraction  of  the 
real  line,  each  has  a  unique  fixed  point  by  the  contraction  mapping  theorem. 
However,  under  the  union  of  these  maps  W  =  wi  u  W2,  the  endpoints  of  the 
intervals  we  constructed,  1/3, 2/3, 1/9, 2/9, 7/9, 8/9, 1/27, . . . ,  also  remain  in 
the  set  after  each  iteration  of  W.  While  we  restricted  our  attention  to  the 
unit  interval,  we  note  that  each  contraction  Wi  is  defined  on  the  real  line  and 
their  union  is  a  contraction  on  the  space  Ji(R).  Hence,  if  we  can  find  a  closed, 
bounded  subset  of  the  real  line  whose  image  under  W  is  itself,  we  know  it  is 
the  unique  fixed  point  of  the  mapping  W.  The  classical  Cantor  set  is  such  a 
point,  because  taking  the  contraction  W  of  the  Cantor  set  again  yields  the 
Cantor  set.  To  see  this,  recall  that  in  our  construction  of  the  Cantor  set  each 
subinterval  contains  a  smaller  copy  of  the  whole  set  because  of  the  infinite 
iteration  we  used.  Likewise,  since  W  maps  the  unit  interval  to  [0,1/3]  and 
[2/3, 1],  each  subinterval  also  contains  a  smaller  copy  of  the  whole  interval. 
Since  the  Cantor  set  is  the  fixed  point  of  this  iterated  function  system  (IFS),  we 
refer  to  it  as  the  attractor  of  W. 

This  illustrates  a  property  of  fractals  known  as  self-similarity.  If  we 
examine  the  interval  [0, 1/3]  in  the  construction  of  the  Cantor  set,  we  see  that 
it  is  a  scaled  copy  of  the  original  unit  interval.  Likewise,  if  we  examine  any 
subinterval  of  the  unit  interval,  we  observe  the  same  result.  This  is  a  unique 
property  of  the  infinite  iteration  scheme  used  to  construct  the  Cantor  set. 
Similarly  it  is  the  result  of  iterating  W  an  infinite  number  of  times.  This 
self-similarity  characteristic  of  fractals  is  very  useful  in  applications  to 
computer  graphics. 


35 


Because  it  was  known  ahead  of  time  that  the  attractor  of  W  would  be 
contained  entirely  within  the  unit  interval,  to  simplify  the  discussion  we 
based  our  previous  analysis  on  mappings  of  the  unit  interval  W:   [0, 1  ]—►[(),  1]. 
Now  consider  the  result  of  mapping  a  different  compact  (i.e.,  closed  and 
bounded)  subset  of  the  real  line.  From  the  Contraction  Mapping  Theorem, 
no  matter  which  nonempty  compact  subset  of  the  real  line  we  begin  with, 
infinite  iteration  of  W  always  produces  the  Cantor  set.  Although  we  proved 
this  principle  formally,  it  certainly  is  not  intuitively  obvious.   In  the  next 
section  we  illustrate  this  idea  with  an  example  in  R2,  where  the  geometry  can 
be  better  appreciated. 

The  question  now  arises,  what  happens  to  the  entire  real  line  under 
iteration  by  W?  It  is  clear  that  W(-°°,  ©©)  =  (-©«,  <*>),  which  appears  to  be  a  fixed 
point  of  W.  However,  the  set  (-<»,  <*>)  is  not  in  the  space  H(R)  since  points 
in  the  space  Ji(R)  are  closed  and  bounded  subsets  of  R,  but  (-~,  <*>)  is  not 
bounded. 

This  thorough  examination  of  the  Cantor  set  as  a  fractal  produced  by 
the  iterated  function  system  W  provides  the  basis  for  similar  results  in 
higher  dimensions.   The  principles  remain  unchanged  and  adding  a 
dimension  yields  ever  more  interesting  results. 

2.      Fractals  in  Two  Dimensions 

We  now  lift  the  setting  from  the  real  line  to  the  Euclidean  plane  R2 
and  the  space  (K(R2),  h(d))  associated  with  it.  Consider  the  following  three 
contractions  of  the  plane: 
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Let  W  =  wi  u  W2  u  W3.  Then  W  is  a  contraction  on  Ji(R2),  and  has  a  unique 
fixed  point,  or  attractor.  To  find  the  attractor,  consider  the  action  of  W  on 
the  unit  right  triangle  T  with  vertices  (0,  0),  (1,  0),  and  (0, 1).  Since  T  is  a 
closed  and  bounded,  nonempty  subset  of  R2,  it  is  a  point  in  H(R2).  Each  of 
the  maps  wi,  W2,  and  W3,  shrinks  any  point  in  R2  towards  the  origin  by  a 
factor  of  1/2  in  the  directions  of  each  coordinate  axis.  The  map  W2  shifts 
each  point  1/2  unit  to  the  right  as  well,  while  W3  shifts  the  result  1/2  unit 
upwards. 

Figure  3.8  shows  T  and  the  result  of  W(T).  Notice  that  the  result  of  W(T)  is 
to  remove  a  "middle  third"  triangle  from  the  original  triangle  T.   Again, 
iterating  T  under  W  yields  the  results  shown  in  Figure  3.9,  displaying  the 
remaining  triangles  from  their  previous  iterations  with  the  middle  thirds 
removed. 

Continuing  to  iterate  W,  we  see  that  the  fixed  point  of  W  is 
obtained  by  continually  removing  the  middle  third  of  each  subtriangle  from 
the  original  triangle  T.  This  attractor  is  shown  in  Figure  3.10  and  is  called  the 
Sierpinski  triangle  (or  Sierpinski  gasket). 


37 


1  x1 

The  Unit  Triangle  T  W(T) 

Figure  3.8  The  effect  of  W  on  the  unit  triangle. 
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Figure  3.9  The  effect  of  W2  on  the  unit  triangle. 

The  same  results  observed  for  the  Cantor  set  in  one  dimension  apply 
to  the  Sierpinski  triangle  in  two  dimensions.   The  Sierpinski  triangle  has 
infinite  detail,  and  is  self  similar:  each  subtriangle  contains  a  scaled  copy  of 
the  original  unit  triangle.   Furthermore,  while  we  used  the  unit  right  triangle 
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T  as  the  starting  point,  by  the  Contraction  Mapping  Theorem  any  point  in  the 
space  H(R2)   converges  to  the  Sierpinski  triangle  under  infinite  iteration  of 
W.   Figure  3.11  shows  how  the  unit  square  with  infinite  iteration  under  W 
also  yields  the  Sierpinski  triangle. 

Figure  3.12  shows  how  the  point  (1, 1)  starts  away  from  the  attractor, 
but  under  iteration  of  W  eventually  converges  to  points  on  the  attractor. 
Figure  3.13  shows  the  effect  of  W6(S)  on  the  unit  square. 


Figure  3.10  The  Sierpinski  triangle. 
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THE  UNIT  SQUARE 
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Figure  3.11  Creating  the  Sierpinski  triangle  by  iterating  the  unit  square. 
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Figure  3.12  The  orbit  of  the  point  (1, 1)  under  iteration  of  W. 


Figure  3.13  The  effect  of  W6(S)  on  the  unit  square. 
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3.      Condensation  Sets 

So  far  we  have  only  considered  examples  of  contractions  where 
0  <  s  <  1.  We  will  now  consider  contractions  where  the  contractivity  factor  is 
s  =  0.   Such  contractions  are  called  condensations,  and  in  one  dimension  this 
corresponds  to    f  (x)  =  c  where  c  is  a  constant.  It  is  clear  that  no  matter  what 
two  points  we  start  with,  the  distance  between  their  first  iterates  under  f  is  0, 
since  c  -  c  =  0,  so  that  f  is  a  contraction. 

In  the  space  H(X),  wo:  H(X)-»H(X)  is  a  condensation  if  wo(x)  =  B  for 
all  x  in  H(X),  where  B  is  some  fixed  point  in  K(X).  If  we  take  a  contraction 
mapping  of  the  space  H,  W  =  wi  u. .  .u  Wn,  and  form  the  union  of  it  with  wo, 
where  wo  is  a  condensation  such  that  W'=  W  u  wo,  then  W  will  be  a 
contraction  with  the  same  contractivity  factor  as  W.  To  see  this,  consider  that 
the  contractivity  factor  of  W  is  s  =  max{si,  S2,  •  • . ,  sn},  and  that  the 
contractivity  factor  of  W  is  s'  =  max{0,  si,  S2, . . . ,  Sn)  =  s. 

Although  the  contractivity  factors  of  the  maps  W  and  W  are  the 
same,  their  dynamics  are  quite  different.  Since  every  iteration  of  W  maps 
points  to  the  point  B,  the  attractor  of  this  iterated  function  system  can  be 
thought  of  as  the  infinite  set  of  every  iteration  of  B  under  W. 

As  an  example,  consider  the  one-dimensional  case  with  the  two 
functions  wi(x)  =  x/2  and  wo(x)  =  1  on  the  real  line.  The  map  wi  is  a 
contraction,  and  wo  is  a  condensation,  and  their  fixed  points  are  0  and  1 
respectively.  If  we  start  with  a  point  on  the  real  line,  for  example  0,  and 
iterate  it  under  W  -  w0  u  wi,  we  see  that  W(0)  =  {0, 1},  W2(0)  ■  {0, 1, 1/2}, 
W3(0)  =  {0, 1, 1/2, 1/4},  and  WHO)  =  {0, 1, 1/2, 1/4, ... ,  l/2n-1}.  Clearly,  this 
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iterated  function  system  generates  an  attractor  which  consists  of  every  iterate 
of  the  condensation  set  {1}  under  the  contraction  wi(x)  =  x/2. 

4.      A  Fractal  Tree 

We  now  use  the  concept  of  condensation  sets  to  construct  a  fractal 
image.  We  first  let  the  set  B  =  {(xi,  X2):  -0.1  <  xi  <  0.1;  0  <  X2  <  1}.  In  other 
words,  B  is  the  filled  rectangle  with  height  1  and  width  .2  with  base  centered 
at  the  origin.   Now  define  the  transformations 


w. 


ft]-- 


w. 


.75cos(tc/4)  -.75sin(7c/4) 
.75sin(7c/4)  .75cos(tc/4) 


EW 


.75cos(-jc/4)  -.75sin(-7c/4) 
.75sin(-Jc/4)  .75cos(-7i/4) 


X2 


Ml 


The  transformations  wi  and  W2  shrink  any  image  by  three-fourths  in  the 
direction  of  both  axes,  wi  rotates  it  45*  clockwise  and  W2  rotates  it  45* 
counterclockwise,  and  then  both  transformations  shift  the  image  up  one  unit. 
The  transformation  wo  is  a  condensation  that  maps  any  point  to  the  set  B. 
If  we  let  W  =  wo  u  wi  u  W2,  and  we  iterate  the  origin  under  W,  we 
see  that  W(0,  0)  yields  the  set  B,  W2  yields  the  set  B  together  with  the  two 
condensed,  rotated,  shifted  copies  of  B,  and  W3  yields  the  set  B  together 
with  two  condensed,  rotated,  and  shifted  copies  of  the  image  of  W2.  If  we 
continue  iterating,  we  see  that  every  iteration  yields  the  set  B  together  with 
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two  condensed,  rotated,  and  shifted  copies  of  the  previous  iteration.  The  first 
seven  iterations  of  the  origin  under  W  are  shown  in  Figure  3.14,  which  is 
the  beginning  of  a  fractal  "tree." 


Figure  3.14  A  fractal  tree. 
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Another  way  of  generating  a  fractal  tree  is  through  an  iterated 
function  system  without  condensation.   If  we  consider  the  transformations 
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,  1  cos(rc/4)  -.25sin(rc/4) 
.lsin(jc/4)    .25cos(tc/4) 


x2 


;H: 


w3 
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.  1  cos(-tu/4)  -.25sin(-Jc/4) 
.  1  sin(-7c/4)    .25cos(-7c/4) 


[:;]*[: 


and  we  consider  the  effect  of  W  =  wi  u  W2  u  W3  on  the  unit  square,  we  see 
the  result  shown  in  Figure  3.15.  If  we  continue  to  iterate  W,  we  produce 
the  attractor  of  this  iterated  function  system,  which  is  a  fractal  tree  with  a 
similar  shape  to  the  one  in  Figure  3.14.  However,  unlike  the  tree  with 
condensation,  the  trunk  and  branches  of  this  tree  are  not  solid,  in  fact,  each 
contains  a  condensed  copy  of  the  entire  attractor.  These  detailed  images  are 
not  only  aesthetically  pleasing,  but  they  lead  to  an  important  application  of 
fractals.  To  store  these  images  in  a  computer  using  traditional  means  would 
require  hundreds  of  lines  of  code;  however,  we  have  used  just  18  numbers 
(three  two-by-two  matrices  and  three  two-vectors)  to  record  all  of  the 
information  required  to  produce  these  images.  This  efficient  way  of  recording 
information  has  led  to  important  applications  of  fractal  geometry,  as  the  next 
section  shows. 
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Figure  3.15  The  action  of  W  on  the  unit  square. 

F.    APPLICATIONS  OF  FRACTALS  TO  COMPUTER  GRAPHICS 

Pictures  that  are  infinitely  detailed  and  self  similar  can  be  stored  in  a 
computer  with  very  small  amounts  of  information.   For  example,  in  the  case 
of  the  fractal  tree  in  the  last  section,  only  18  numbers  were  required.  But  what 
if  an  image  is  not  self  similar?  Is  it  still  possible  to  use  contraction  mappings 
to  store  the  information?  The  remarkable  answer  is  "yes."  Perhaps  even 
more  surprising,  this  procedure  can  be  achieved  to  any  desired  degree  of 
accuracy. 

A  formal  description  of  using  contractive  iterated  function  systems  to 
store  graphical  information  can  be  found  in  Barnsley  (1988)  or  Falconer  (1990). 
Instead  of  presenting  the  proofs  here,  we  give  an  intuitive  treatment  of  the 
technique  through  an  example.  Then  we  show  how  any  desired  degree  of 
accuracy  can  be  attained. 
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Recall  that  affine  transformations  can  have  many  different  actions  on  a 
set,  including  shearing,  shrinking  or  stretching,  rotating,  and  shifting. 
Suppose  we  have  an  image  we  wish  to  reproduce;  let  that  be  our  attractor. 
We  then  search  for  affine  transformations  that  could  paste  distorted  copies  of 
the  set  back  on  top  of  itself  to  approximate  the  original  image  as  closely  as 
possible.  This  is  done  through  the  "collage  theorem"  due  to  Barnsley.   We 
present  an  example  here. 

Consider  the  image  in  Figure  3.16.  While  this  figure  is  not  self  similar, 
we  can  approximate  it  by  making  smaller  copies  of  it  (through  contractions) 
and  "pasting"  them  back  on  itself  (through  affine  transformations)  to  cover, 


Figure  3.16  A  nonself-similar  image. 
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as  much  as  possible,  the  original  figure.  One  such  covering  is  shown  in 
Figure  3.17.   The  four  corners  of  the  original  image  and  their  movement 
under  the  coverings  are  indicated  by  the  letters  A  through  D.  Note  that  a 
certain  amount  of  overlapping  of  the  covering  is  required,  for  if  this  were 
just-touching  or  totally  disconnected,  the  attractor  would  be  too  sparse  to 
represent  the  original  image. 
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Figure  3.17  Covering  an  image  with  transformed  copies  of  itself. 


Recalling  that  an  affine  transformation  is  uniquely  determined  by  its 
action  on  any  three  noncollinear  points,  we  could  easily  determine  the  actual 
transformations  required  to  produce  the  seven  mappings  indicated  in  the 
figure.  Under  repeated  iteration  of  this  iterated  function  system,  we  could 
then  reproduce  the  image  in  a  similar  way  as  that  used  to  produce  other 
attractors,  with  only  42  pieces  of  data,  since  each  affine  transformation 
requires  six  numbers. 
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But  notice  that  the  attractor  of  this  iterated  function  system  is  not  the 
original  image,  since  its  boundaries  are  infinitely  detailed,  as  opposed  to  the 
solid  boundaries  in  the  original  image  (see  Figure  3.18).   This  "fuzziness" 
seems  to  contradict  the  earlier  claim  that  we  could  approximate  a  computer 
image  to  any  degree  of  accuracy.  But  the  above  system  used  only  seven 
transformations.  Noting  that  any  computer  image  is  displayed  by  a  finite 
number  of  pixels,  which  are  either  on  or  off,  we  could  push  this  technique  to 
its  extreme  to  reproduce  exactly  the  original  image.  That  is,  we  could 
condense  the  entire  set  down  to  the  size  of  one  pixel,  and  then  shift  that  pixel 
through  affine  transformations  to  every  pixel  in  the  original  image.   This  will 
indeed  give  an  exact  reproduction  of  the  original  image.  However,  it  would 
require  six  numbers  for  every  pixel  in  the  image. 


Figure  3.18  The  attractor  with  a  fractal-like  boundary. 
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While  this  technique  obviously  defeats  the  concept  of  efficient  storage  of 
information,  it  nevertheless  shows  that  we  can  achieve  any  desired  degree  of 
accuracy.   Depending  on  the  requirements  of  a  given  project,  iterated  function 
systems  and  their  attractors  can  frequently  be  used  to  achieve  the  desired 
degree  of  accuracy  much  more  efficiently  than  the  traditional  methods  of 
storing  information. 

So  far  we  have  claimed  that  only  the  information  contained  in  the  affine 
transformations  is  needed  to  use  this  technique.  In  fact,  we  also  require  a 
program  to  produce  the  image  based  on  the  input.  However,  these  are  quite 
common  and  there  are  a  variety  of  algorithms  that  produce  the  desired 
results.  Thus,  in  a  situation  where  many  different  figures  must  be  recorded,  a 
single  program  will  suffice  for  creating  the  graphics  once  the  transformations 
are  stored.  It  follows  that  this  method  is  more  efficient  than  many  traditional 
computer  graphics  techniques. 

G    THE  ADDRESSES  OF  POINTS  ON  FRACTALS 

A  useful  technique  in  analyzing  fractals  is  to  address  every  point  on  an 
attractor  by  the  sequence  of  transformations  that  led  to  that  point  being  in  the 
attractor.  This  leads  to  the  classification  of  iterated  function  systems  as  being 
"totally  disconnected,"  "just  touching,"  or  "overlapping,"  and  helps  us 
analyze  chaotic  dynamics  on  fractals  in  Chapter  IV. 

We  will  start  with  transformations  of  the  real  line  to  illustrate  this 
concept  in  the  simplest  setting  possible.  Consider  the  system  wo(x)  =  x/3  and 
W2(x)  =  x/3  +  2/3  (the  reason  for  selecting  these  subscripts  will  soon  become 
apparent).  We  already  know  that  repeated  iteration  of  W  =  wquw2  yields 


50 


the  Cantor  middle-thirds  set  as  an  attractor.   (Recall  that  there  are  no  intervals 
in  the  Cantor  set,  and  that  the  intersection  of  wo  and  W2  applied  to  the  unit 
interval,  or  any  subset  thereof,  is  empty).   Considering  the  iteration  of  W  on 
the  unit  interval,  we  now  begin  to  address  the  points  in  its  attractor. 

Under  the  first  iteration  of  W,  the  unit  interval  is  mapped  to  the  interval 
[0,1/3]  by  wo  and  to  [2/3,1]  by  W2-  Hence,  we  begin  the  address  of  every 
point  in  [0, 1/3]  with  0,  and  we  begin  the  address  of  every  point  in  [2/3, 1] 
with  2.  Under  the  second  iteration  of  W,  the  interval  [0,1/3]  is  mapped  to 
[0,1/9]  by  wo  and  to  [2/3,7/9]  by  W2,  whereas  the  interval  [2/3,1]  is 
mapped  to  [2/9,1/3]  by  wo  and  to  [8/9,1]  by  W2-  Hence,  the  second  number 
in  the  address  of  all  points  in  [0,1/9]  and  [2/3,7/9]  is  0,  and  the  second 
number  in  the  address  of  all  points  in  [2/9, 1/3]  and  [8/9, 1]  is  2.  Figure  3.19 
shows  the  first  three  steps  of  this  process,  and  the  beginning  of  the  resulting 
addresses  of  every  point  in  the  Cantor  set. 
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Figure  3.19  Addressing  the  attractor  of  an  iterated  function  system. 

It  is  clear  that  under  infinite  iteration,  every  point  in  the  Cantor  set  has  a 
unique  address  consisting  of  an  infinite  string  of  0s  and  2s.  Note  that  the 
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origin  has  the  address  000. . . ,  and  that  1  has  the  address  222. . . .  What  we 
have  by  the  convenient  selection  of  0  and  2  as  the  subscripts  is  that  every 
point  in  the  Cantor  set  is  addressed  by  its  ternary  expansion,  as  described  in 
Chapter  n. 

The  addressing  scheme  used  here  is  not  unique,  as  we  could  have  also 
chosen  0  and  1  as  subscripts  resulting  in  every  point  in  the  Cantor  set 
having  an  address  from  the  code  space  £2-  Again,  note  that  every  point  in 
£2  would  be  an  address  of  a  point  in  the  Cantor  set,  and  that  any  dynamics  we 
performed  on  the  code  space  X2  could  easily  be  applied  to  the  Cantor  set. 
This  one-to-one  correspondence  between  code  space  and  fractals  is  very  useful 
in  the  discussion  of  chaotic  dynamics  on  attractors  of  iterated  function 
systems.  (As  an  aside,  we  could  have  just  as  easily  selected  1  and  2  as 
subscripts  to  demonstrate  further  that  an  addressing  scheme  is  not  unique. 
Nevertheless,  all  of  these  systems  have  a  one-to-one  correspondence  between 
them.) 

Notice  that  every  point  in  the  Cantor  set  has  a  unique  address  associated 
with  it.  This  characteristic  allows  us  to  classify  the  iterated  function  system 
W  =  wo  u  W2  as  being  totally  disconnected.  While  the  attracting  set  (the 
Cantor  set)  is  disconnected,  we  associate  the  classification  "totally 
disconnected"  not  with  the  attractor,  but  rather  with  the  iterated  function 
system  that  produced  it.  The  reason  for  this  will  become  apparent  when  we 
look  at  "just  touching"  and  "overlapping"  iterated  function  systems  which 
have  the  same  set  as  their  attractor. 

Now  consider  the  iterated  function  system  W  =  wi  u  W2,  where 
wi(x)  =  x/2  and  W2(x)  =  x/2  +  1/2.  The  unique  fixed  point  (and  hence  the 
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attractor)  of  this  system  is  clearly  the  unit  interval.  So  if  we  consider  the 
action  of  this  system  on  [0, 1],  we  can  build  addresses  for  every  point  in  that 
interval.  Starting  with  [0 ,1],  we  see  that  wi  maps  it  to  [0,1/2]  so  we  assign 
to  that  interval  the  first  address  of  1,  and  similarly  points  in  [1/2,1]  have  the 
first  address  2.  Already  we  can  see  that  there  is  some  ambiguity  with  the 
point  1/2,  since  it  appears  to  be  receiving  two  distinct  addresses.  Continuing, 
we  see  that  points  in  the  interval  [0,1/4]  have  addresses  starting  with  11, 
while  points  in  the  intervals  [1/4,1/2],  [1/2,3/4],  and  [3/4,1]  have  addresses 
12,  21,  and  22,  respectively.  Now  the  points  1/4  and  3/4  have  been  added 
to  the  list  of  points  with  dual  representations.  Figure  3.20  shows  the  first 
three  iterations  of  this  addressing  scheme. 

0  1/4  1/2  3/4  1 


11 HO  21 

12        — — —        22 
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Figure  3.20  Addressing  the  attractor  of  a  just-touching  IFS. 

Again,  it  is  clear  that  every  point  in  the  unit  interval  receives  an  address, 
but  points  of  the  form  i/2n,  i  =  1, 2, 3, . . . ,  2n  -  1  have  dual  representations 
under  this  scheme.  The  set  of  points  where  the  addresses  touch,  although 
infinite,  is  countable  (i.e.,  is  in  one-to-one  correspondence  with  the  natural 
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numbers)  since  at  every  step  we  add  a  finite  number  of  points  with  dual 
addresses.  Also,  notice  that  no  intervals  have  dual  addresses,  since  between 
every  pair  of  points  with  dual  addresses,  there  are  points  with  unique 
addresses.  Iterated  function  systems  with  this  characteristic  are  classified  as 
just-touching.  This  fits  the  intuitive  concept  of  this  definition  since  the 
results  of  iterating  the  attractor  [0, 1]  under  wi  and  W2  results  in  two  intervals 
which  "just  touch"  at  one  point. 

Finally,  we  consider  the  system  W  =  wjUW2  where  wi(x)  =  2x/3  and 
W2(x)  =  2x/3  +  1/3.  Again,  the  attractor  of  this  system  is  clearly  [0, 1],  but 
when  we  address  points  in  the  attractor,  we  encounter  considerable 
ambiguity.  Under  the  first  iteration  of  W  applied  to  [0, 1],  we  see  that  the 
interval  [0,2/3]  gets  the  address  1,  while  the  interval  [1/3,1]  gets  the  address 
2.  Here  an  entire  interval  [1/3,2/3]  of  points  has  a  dual  address  after  just  the 
first  iteration.  Moreover,  at  every  step  we  add  intervals  with  a  similar 
ambiguity.  The  first  three  steps  in  this  addressing  scheme  are  shown  in 
Figure  3.21.  Notice  that  in  order  to  remove  the  ambiguity  of  the  overlapping 
addresses,  we  have  to  "lift"  the  unit  interval  into  a  second  dimension  as  the 
figure  shows. 

The  ambiguity  in  this  addressing  system  is  much  greater  than  for  the  just- 
touching  iterated  function  system.  In  fact,  at  every  step  an  uncountable 
number  of  points  with  multiple  addresses  is  obtained.  Here,  the  number  of 
points  with  multiple  addresses  far  outnumbers  the  points  with  unique 
addresses.  What  is  more,  for  this  system,  the  only  points  with  unique 
addresses  are  0  and  1.  This  feature  leads  to  the  intuitive  definition  of  an 
overlapping  iterated  function  system.  Notice  that  each  iteration  of  wi  and 
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W2  results  in  an  interval  where  different  addresses  overlap.   Because  the 
attractor  of  this  system  is  the  same  as  in  the  example  of  a  just-touching 
iterated  function  system,  it  should  now  be  clear  why  the  classifications 
"totally  disconnected,"  "just  touching,"  and  "overlapping"  are  applied  to  the 
iterated  function  system  rather  then  to  the  attractor  itself. 
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Figure  3.21  Addressing  the  attractor  of  an  overlapping  IFS. 

With  the  concept  of  addressing  points  on  attractors  and  classifying 
iterated  function  systems  in  a  simple  one-dimensional  setting,  let  us  now 
turn  to  a  two-  dimensional  example.   Recall  the  iterated  function  system  that 
produced  the  Sierpinski  triangle  as  an  attractor: 
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where  W  =  wi  u  W2  u  W3.  Iterating  the  unit  right  triangle  under  W  yields  the 
Sierpinski  triangle  by  removing  the  "middle-thirds  triangle"  at  every  step. 
Looking  at  just  one  iteration  of  W,  we  see  that  this  is  a  just-touching  iterated 
function  system:  the  image  of  wi  intersects  the  image  of  W2  at  (1/2,  0),  the 
images  of  W2  and  W3  intersect  at  (1/2, 1/2),  and  the  images  of  wi  and  W3 
intersect  at  (0, 1/2).  An  addressing  scheme  for  the  Sierpinski  triangle  is 
shown  in  Figure  3.22. 


Figure  3.22    An  addressing  scheme  for  the  Sierpinski  triangle. 

As  before,  the  vertices  created  at  each  step  have  two  distinct 
addresses.  For  example,  the  point  (1/2, 1/2)  has  the  addresses  2333. . .  and 
3222.  .  .  .  The  exceptions  are  the  three  original  vertices  which  have  the  unique 
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addresses  111. . . ,  222. . . ,  and    333. . . .  Again,  we  could  have  chosen  the 
addresses  to  consist  of  Os,  Is,  and  2s,  to  create  a  correspondence  between  the 
Sierpinski  triangle  and  code  space  X3.  However,  it  is  just  as  convenient  to 
create  a  modified  code  space  X3  to  be  {X1X2X3. . . :  xj  e {1, 2, 3}}.  The  metric 
associated  with  this  space  £'3  is 

d(x,y)  =  £    !    /', 
£       41 

so  Z'3  is  clearly  in  one-to-one  correspondence  with  the  Sierpinski  triangle. 

A  more  difficult  addressing  problem  occurs  with  an  overlapping 

iterated  function  system  in  two  dimensions.   Recall  that  in  an  overlapping 

iterated  function  system  of  the  real  line,  we  had  to  "lift"  the  attractor  into  a 

second  dimension  to  uniquely  determine  addresses.  Similarly,  in  two 

dimensions,  we  must  lift  the  attractor  into  a  third  dimension  to  address  its 

points  uniquely.  For  example,  consider  the  overlapping  iterated  function 

system  of  the  plane  W  =  wi  u  W2  u  W3,  where 


w 


The  attractor  for  this  iterated  function  system  is  the  "filled" 
Sierpinski  triangle,  or  simply  the  unit  right  triangle.   Figure  3.23  shows  the 
first  iteration  of  W  and  how  the  iterates  of  wi,  W2,  and  W3  overlap.  Any 


57 


point  along  the  xi -axis  in  the  interval  [1/3,2/3]  has  an  address  beginning 
with  a  1  or  a  2,  and  entire  regions  of  the  plane  have  similar  ambiguities. 


1  or  3 


3V 

^v2or3 

1 

C\ 

1  or  2 
Figure  3.23  Ambiguous  addresses  for  the  attractor  of  an  overlapping  IFS. 

While  points  on  the  attractor  always  correspond  to  multiple 
addresses,  we  would  like  to  identify  every  possible  address  with  a  unique 
point.  This  is  useful,  for  example,  when  considering  the  dynamics  of  the  shift 
map  on  code  space  associated  with  the  attractor  of  an  iterated  function  system. 
Developing  such  an  addressing  scheme  in  two  dimensions  for  this  attractor 
would  be  extremely  difficult.  However,  "lifting"  the  attractor  into  a  third 
dimension  allows  us  to  create  an  unambiguous  address  system  for  the 
attractor.  Figure  3.24  shows  how  the  first  iterate  of  W  (from  Figure  3.23)  is 
lifted  into  the  third  dimension.   This  removes  the  ambiguous  regions  from 
Figure  3.23. 

We  now  lift  the  entire  attractor  into  the  third  dimension.   To 
accomplish  this,  cross  the  Euclidean  plane  with  the  unit  interval  [0, 1]  which 
is  represented  by  its  base  four  expansion  so  that  every  point  in  the  interval  is 
represented  by  a  number  .X1X2X3. . .  where  Xi  €  {0, 1, 2, 3}.  By  associating  the 
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address  of  each  iterate  of  W  with  its  base  four  expansion,  we  see  that  we  get  a 
Cantor-type  set  construction  along  the  vertical  axis  when  viewed  from  the 
side.  Since  we  used  the  subscripts  1,  2,  and  3  for  the  iterated  function  system, 
every  number  in  [0, 1]  that  does  not  have  a  0  in  its  base  four  expansion  will 
have  an  iterate  associated  with  it.  The  result  of  this  process,  taken  to  three 
iterates,  is  shown  in  figure  3.25.  Having  an  addressing  scheme  such  as  this  to 
locate  a  unique  point  proves  quite  useful  when  we  study  chaotic  dynamical 
systems  on  fractals  in  Chapter  IV. 


Figure  3.24  Lifting  an  attractor  into  the  third  dimension. 
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Figure  3.25  Addressing  the  lifted  attractor  of  an  overlapping  DFS. 

H.    FRACTAL  DIMENSION 

One  of  the  most  useful  concepts  in  the  application  of  fractal  geometry  is 
that  of  fractal  dimension.  Fractal  dimension  provides  a  measure  of  the  size 
or  "dimension"  of  an  object,  whether  it  is  the  attractor  of  an  iterated  function 
system,  a  more  familiar  geometric  figure,  or  something  arrived  at  through 
collection  of  real-world  data.  Fractal  dimension  has  worked  its  way  into 
many  fields  such  as  physics,  meteorology,  aeronautical  engineering,  and 
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oceanography,  and  we  provide  some  examples  of  these  applications  shortly. 
To  begin  we  develop  this  concept  mathematically. 

The  fractal  dimension  of  an  object  is  the  assignment  of  a  number  to  the 
object  that  represents  how  much  space  the  object  takes  up  in  its  ambient 
space.  We  are  all  familiar  with  objects  in  one,  two,  or  three  dimensions.  That 
an  object  may  have  a  non-integer  or  fractional  dimension  seems 
counterintuitive.   Nevertheless,  you  will  see  that  the  definition  of  fractal 
dimension  supports  completely  our  intuitive  concept  of  one,  two  or  three 
dimensional  objects.   With  this  in  mind,  we  give  some  preliminary 
definitions  that  are  needed  to  develop  the  concept  of  fractal  dimension. 

An  e-ball  about  a  point  xo  in  a  metric  space  is  the  set  of  all  points  in  the 
space  within  distance  e  >  0  from  the  point  xo-  Conventional  notation  for  an 
e-ball  is  B(xo,  e)  =  {x:  d(x,  xo)  <  e}.  Notice  that  by  using  the  "<"  sign  in  the 
definition,  the  e-balls  in  our  discussion  are  closed.  While  this  notation  may 
be  slightly  unfamiliar,  the  concept  certainly  is  not.  For  if  we  take  the  real  line 
with  the  standard  distance  function  as  our  metric  space,  then  B(l,  1/2)  is 
simply  the  closed  interval  [0.5, 1.5].  Furthermore,  if  the  metric  space  is  the 
Euclidean  plane  R2,  then  B((0,  0),1)  is  simply  the  closed  unit  disc  centered  at 
the  origin. 

Given  a  closed  bounded  nonempty  subset  of  a  metric  space,  we  want  to 
cover  that  subset  with  balls  of  specified  radius  e.  In  the  case  of  the  real  line,  to 
cover  the  unit  interval  with  balls  of  radius  1/4  requires  at  least  two  (centered 
at  1/4  and  3/4).  We  could  use  any  number  of  such  balls  if  we  allow 
overlapping.  We  are  primarily  interested  in  covering  a  set  with  the  smallest 
number  of  balls  possible,  so  we  define  the  integer  N(A,  e)  to  be  the  smallest 
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number  of  closed  balls  of  radius  e  needed  to  cover  the  set  A.  Specifically,  for 

any  set  A  in  our  metric  space,  we  define  N(A,  e)  to  be  the  smallest  positive 

integer  M  such  that 

M 

u   B(xn,  e)z>A, 
n=l 

for  some  distinct  set  of  points  {xn:  n  =  1, 2, . . . ,  M}  in  the  metric  space. 

The  Heine-Borel  Theorem  (see  an  advanced  calculus  text)  guarantees  that 
we  are  always  able  to  cover  a  closed  bounded  set  with  a  finite  number  of 
e-balls.  Since  every  set  of  positive  integers  has  a  smallest  member,  the 
number  N(A,  e)  uniquely  exists. 

We  are  now  prepared  to  define  fractal  dimension.  If  we  let  A  be  a  point 

in  the  space  H(X)  (i.e.,  A  is  a  closed  bounded  nonempty  subset  of  the  metric 

space  X)  the  quantity  D  is  defined  as 

D  =  limln(N(A,e)) 
e-*0    ln(l/e) 

The  number  D,  if  it  exists,  is  called  the  fractal  dimension  of  A  and  is  denoted 

D(A). 

Let  us  clarify  this  concept  with  some  intuitive  examples.  First,  consider  a 
point  x  in  R2.  No  matter  how  small  we  choose  e,  we  can  always  cover  the 
point  x  with  a  single  ball  of  radius  e,  so  that  D(x)  =  lime^0ln(N(x/  e))/ln(l/e) 
=  lime-xjlndVlnd/e)  =  0.  Hence,  for  any  point  x  in  R2,  we  see  that  its  fractal 
dimension  is  zero  (which  fits  our  intuitive  concept  of  a  point  having  zero 
dimension). 

Consider  next  the  unit  interval  [0, 1]  as  a  subset  of  the  real  line. 
Recalling  that  in  R,  e-balls  are  the  closed  intervals  of  length  2e,  we  begin  to 
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cover  [0, 1]  with  e-balls.  If  we  let  e  =  1/2, 1/4, 1/8, .. .  ,l/2n, . . .  then  the 
number  of  intervals  required  to  cover  [0, 1]  are  1,  2, 4, ... ,  and  2n"1,  respec- 
tively.   Clearly,  n-»°°  implies  e->0,  so  D([0, 1])  =  lime^oln(N(A,  e))/ln(l /e)  = 
limn_>0oln(2n-1)/ln(2n)  =  limn_>0o(n-l)ln2/nln2  =  limn^oo(n  -  l)/n  =  1.  Hence, 
the  fractal  dimension  of  a  line  segment  is   1,  which  is  consistent  with  our 
ordinary  concept  of  dimension. 

In  the  previous  example  we  were  very  specific  in  the  way  we  let  e-»0. 
This  example  leads  us  to  an  equivalent  definition  of  fractal  dimension.   In  the 
same  setting  as  before,  let  en  =  rn  for  0  <  r  <  1  and  n  =  1, 2, 3, . . .  Then 

D(A)  =  limn_^ln(N(A,  en))/ln(l/en). 
To  show  that  these  definitions  are  equivalent,  let  f(e)  =  Max{en:  en  <  e},  and 
assume  that  e  <  r.  Then  f(e)  £  e  <  f(e)/r,  and  N(A,  f(e»  >  N(A,  e)  £ 

N(A,  f(e)/r).  Since,  for  x  >  1,  ln(x)  is  a  positive  increasing  function,  we  have 

ln(N(A,  f(e)/r))  ^  ln(N(A,  e))  ^  ln(N(A,  f(e))) 

ln(l/f(e))  ln(l/e)  ln(r/f(e)) 

We  can  assume  that  N(A,  £)-»«>  as  e-»0,  or  the  result  follows  immediately. 
Considering  the  left  side,  and  taking  the  limit  as  e-»0,  since  for  some  n, 
f (e)  =  en,  we  have 

lim["ln(N(A,  f(e)/r))1     Km  |to(N(A,  Bn_j)) 


£->°L     ln(l/f(e)) 


lnfl/tej 


liml"    1"(N(A,  £„,!»    1     Bm|"lii(N(A,e,ft' 


H" 


lnd/ej 


Similarly,  the  right  side  yields 


lim 

£->C 


ln(N(A,f(e)))"j_  KmT'nWA.eJ) 
ln(r/f(e))        =n>\     InGAg 
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lim 


ln(N(A>£n)) 
ln(r)+  ln(l/en) 


lim 

n->« 


ln(N(A,en)) 
lnd/e^ 


The  pinching  theorem  from  calculus  establishes  the  equivalence  of  these 
definitions. 

Notice  that  in  R2  with  the  Euclidean  metric,  the  e-balls  always  give 
closed  discs.  Frequently  we  find  it  more  convenient  to  use  "boxes"  or  closed 
squares  to  cover  our  sets.  For  example,  if  we  want  to  cover  the  unit  square 
[0,  1]X[0,  1],  it  would  be  much  easier  to  calculate  the  required  number  of  boxes 
of  a  given  size  than  to  calculate  the  number  of  e-balls.  It  turns  out  that  we  can 
obtain  the  same  value  of  fractal  dimension  using  closed  boxes  as  we  obtain 
using  e-balls. 

First,  consider  a  grid  in  R2.  Figure  3.26  shows  a  covering  of  this  grid  by 
closed  discs,  and  we  see  that  if  we  were  to  refine  our  grid  and  let  the  radii  of 
the  discs  go  to  zero,  then  the  number  of  boxes  and  the  number  of  discs  will 
remain  in  one-to-one  correspondence. 


Figure  3.26  Covering  a  grid  with  closed  discs. 
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The  result  of  placing  finer  and  finer  grids  with  side  length  l/2n  over  a  set 
in  R2,  and  letting  Nn(A)  be  the  number  of  grid  squares  that  intersect  the 
object,  gives  a  third  equivalent  definition  of  fractal  dimension  as  follows: 


D(A)  =  fr" 


ln(N 


n(A))1 
2")| 


ln(2") 

Actually  several  equivalent  definitions  of  fractal  dimension  are  used  in 
different  applications,  depending  on  which  is  the  most  convenient.   All  of  the 
definitions  give  the  same  result.  Therefore,  both  theoretical  and 
experimental  applications  are  quite  easy  to  perform  in  many  different 
situations.  Different  examples  of  these  definitions  for  Euclidean  spaces  are 


ixA)=H"ln(N(A>e))l 

e-°L    ln(l/e) 


where  N(A,  e)  is 

1.  The  smallest  number  of  closed  balls  of  radius  e  that  cover  A  (our  first 
definition); 

2.  The  smallest  number  of  boxes  of  side  e  that  cover  A; 

3.  The  smallest  number  of  sets  of  diameter  of  at  most  e  that  cover  A 
(where  the  diameter  of  a  general  set  is  the  largest  distance  between  all 
pairs  of  points  in  the  set);  and 

4.  The  largest  number  of  disjoint  balls  of  radius  e  with  centers  in  A. 

We  have  presented  these  definitions  in  the  context  of  R  and  R2,  but  they 
are  easily  extended  to  R3,  where  e-balls  become  closed  spheres  and  e-boxes 
become  cubes.  More  generally,  the  same  theory  holds  for  Rn  for  any  integer 
n,  although  the  applications  are  of  a  more  theoretical  nature. 
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Now  that  we  have  discussed  equivalent  definitions  of  fractal  dimension, 
let  us  use  the  most  convenient  one  to  compute  the  dimension  of  the  closed 
unit  square  [0, 1]X[0, 1].  Using  squares  with  side  length  1, 1/2, 1/4, .. .  l/2n  in 
R2,  we  see  that  1, 4, 16, . . . ,  4n  of  them  are  required  to  cover  the  unit  square, 
respectively.   Then 

D(to,  i]x[o,  id = ^n^ 


we  have 


n(A))1 
*    J 


rin(4°)"[ 
Iln(2")J 


D([0,l]Xt0,l])=    fa|iS(£|=    fa^  =  2. 


n^°°|_ln(2^J     n_>0° 
Hence,  our  definition  of  fractal  dimension  again  supports  our  knowledge  that 
the  closed  unit  square  has  dimension  2. 

We  now  apply  the  definition  of  fractal  dimension  to  the  Cantor  set  and 
the  Sierpinski  triangle.   Intuitively,  one  might  believe  that,  since  these  objects 
lie  in  R  and  R2  respectively,  they  should  have  dimensions  1  and  2.  But  recall 
that  we  removed  a  certain  amount  of  length  or  area  at  each  iteration  for  an 
infinite  number  of  iterations.  Thus  it  might  be  plausible  for  these  objects  to 
have  a  smaller  dimension  than  that  of  the  spaces  in  which  they  exist. 

Consider  first  the  Cantor  set  as  a  subset  of  the  unit  interval.  Let  the  e-balls 
in  R  be  the  intervals  of  length  1, 1/3, 1/9, ... ,  l/3n.  One  can  easily  see  that 
1, 2, 4, . . . ,  2n  balls  are  required  to  cover  the  set,  respectively.  Hence,  the 
fractal  dimension  of  the  Cantor  set  is 


lim 
D(Q  = 

n— »<*> 


ln(2") 


ln(3") 


lim    nln2_ln2 
n_»©o  nln3  ~  ln3' 
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so  that  the  Cantor  set  has  dimension  of  approximately  .6309.  We  removed  a 
total  length  of  1   from  the  unit  interval  when  constructing  the  Cantor  set,  but 
still  left  an  uncountable  number  of  points.   Thus  intuitively  we  might  expect 
its  dimension  to  be  somewhere  between  0  and  1.  Fractal  dimension 
provides  a  measure  of  how  many  points  actually  remain  in  the  set. 

Moving  next  to  the  Sierpinski  triangle,  if  we  use  squares  of  side  length    1, 
1/2, 1/4, ... ,  l/2n,  then  1, 3, 9, . . . ,  3n  squares  are  required  to  cover  the  set 
(see  Figure  3.27). 


Length  =  1 


Length  =  1/2 


Length  =  1/4 


One  Square  Three  Squares  Nine  Squares 

Figure  3.27  Covering  the  Sierpinski  triangle  with  grids  of  size  l/2n. 

Hence,  the  dimension  of  the  Sierpinski  triangle  is  given  by 


D(s)  = 


lim 


n— »«> 


ln(3") 


ln(2") 


Km    nln3=ln3 
n_>oo  nln2     ln2' 


and  the  Sierpinski  triangle  has  dimension  of  approximately  1.585.  Again, 
considering  the  amount  of  area  removed  from  the  unit  triangle,  it  seems 
reasonable  that  the  dimension  should  be  between  1  and  2.  It  now  becomes 
clear  why  the  term  "fractal"  was  chosen  to  describe  objects  like  the  Cantor  set 
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and  the  Sierpinski  triangle,  as  they  have  dimensions  which  are  fractions  of 
whole  numbers. 

While  we  have  shown  how  to  determine  the  fractal  dimension 
theoretically  in  some  very  simple  cases,  we  are  not  always  able  to  make  this 
computation  by  just  looking  at  the  attractor  geometrically.  However,  if  the 
iterated  function  system  that  created  the  attractor  is  known,  the  fractal 
dimension  of  the  attractor  can  be  determined  by  analyzing  the  iterated 
function  system  itself.  We  discuss  this  idea  next,  and  state  a  general  result. 

If  we  know  that  W  =  wi  u  W2  u. .  .u  wn  is  the  iterated  function  system 
producing  an  attractor,  and  if  each  wj  is  a  similitude  with  contractivity  factor 
si  (as  defined  in  Section  m.3,  in  that  it  is  a  shrinking  in  each  direction  by  the 
factor  r  such  that  Si  =  n,  followed  by  a  rotation  and  a  translation)  then  if  the 
iterated  function  system  is  totally  disconnected  or  just  touching,  the  attractor 
has  fractal  dimension  D(A)  given  by  the  unique  solution  of 

i=l 

where  D(A)  will  be  between  0  and  the  dimension  of  the  ambient  Euclidean 
space.  If  the  iterated  function  system  is  overlapping,  then  we  only  get  a 
bound  for  D(A)  given  by  D(A)  <  D,  where 

iiSiiD=i. 

i=l 

Even  if  we  do  not  know  the  exact  iterated  function  system  that  produced 
the  attractor,  we  can  measure  the  contraction  of  our  set  under  each 
similitude,  and  solve  for  a  rough  approximation  of  D(A).  While  the  proof  of 
this  result  is  quite  complicated,  (see  Barnsley,  1988,  p.  185)  it  does  give  a  very 
useful  theoretical  tool  for  determining  the  fractal  dimension  of  an  attractor. 
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To  demonstrate  its  use  on  the  Sierpinski  triangle,  which  we  know  is  created 
with  three  similitudes  with  contractivity  factor  of  1/2,  and  knowing  that  the 
Sierpinski  triangle  is  just  touching,  we  have  (1/2)D  +  (1/2)D  +  (1/2)D  =  1, 
from  which  we  get  D  =  In3/ln2,  as  expected.  Unfortunately,  the  theoretical 
determination  of  fractal  dimension  is  not  always  possible,  so  we  must 
frequently  use  experimental  methods  in  applications. 

L     EXPERIMENTAL  DETERMINATION  OF  FRACTAL  DIMENSION 

Fractal  dimension,  as  previously  discussed,  is  perhaps  the  most  useful 
aspect  of  fractal  geometry  in  real-world  applications.  It  has  permitted  the 
analysis  of  numerous  natural  phenomena  that  previously  were  inaccessible 
to  scientists  because  of  a  lack  of  tools  and  theory.  In  nature,  nice  standard 
geometric  shapes  are  primarily  the  exception  rather  than  the  rule;  irregular, 
fractal-like  shapes  are  found  everywhere.  Consider,  for  instance,  the  shapes 
of  clouds,  coastlines,  geological  formations,  atmospheric  phenomena,  plants, 
crystals,  and  microorganisms.  None  of  these  objects  can  be  described  by  a 
"nice"  geometric  shape,  but  fractal  geometry  has  provided  a  new  way  of 
analyzing  them. 

One  way  scientists  have  studied  nature  through  fractal  geometry  is  by 
comparing  the  fractal  dimension  of  an  observed  natural  phenomenon  to  a 
system  with  a  known  fractal  dimension  (whether  the  known  system  was 
created  in  a  laboratory,  created  mathematically,  or  is  a  natural  phenomenon 
that  has  previously  been  studied).  Similar  fractal  dimensions  are  used  to 
support  hypotheses  about  the  similarity  of  the  two  systems.  Thus,  the 
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experimental  determination  of  fractal  dimension  has  become  quite 
important. 

Currently,  there  is  no  definite  theory  for  determining  fractal  dimension 
experimentally.  However,  what  follows  is  one  of  the  most  widely  used 
methods.  When  real-world  data  are  gathered  experimentally,  they  can 
usually  be  plotted  as  a  subset  of  one  of  the  Euclidean  spaces  R,  R2,  or  R3. 
Using  e-balls  of  varying  radii,  one  can  "cover"  the  data  with  the  smallest 
number  of  balls  possible,  obtaining  values  for  e  and  N(A,  e).  This  procedure 
is  not  as  simple  as  it  may  at  first  seem,  because  though  we  know  when  a  set  of 
points  is  covered,  we  can  never  be  sure  that  the  smallest  number  of  balls 
possible  has  been  used.  Nevertheless,  with  the  aid  of  computers,  a  very  good 
approximation  of  N(A,  e)  for  a  given  e  can  frequently  be  obtained. 

In  R2,  one  way  of  accomplishing  this  is  to  superimpose  a  mesh  of  grid 
size  e  over  the  data  points,  and  count  the  number  of  squares  that  contain 
data  points;  this  procedure  has  an  obvious  analogy  in  R3.  By  varying  the 
value  of  e,  pairs  of  values  for  e  and  N(A,  e)  can  be  obtained.  Plotting  the 
logs  of  these  values  against  each  other  over  a  range  of  different  values  of  £ 
should  yield  (approximately)  a  straight  line  from  which  the  fractal  dimension 
of  the  set  can  be  determined. 

There  are  some  drawbacks  to  this  method.  The  first  problem  is  that  the 
theory  we  have  developed  for  fractal  dimension  applies  only  to  closed, 
bounded  subsets  of  some  Euclidean  space.  The  real  world,  on  the  other  hand, 
does  not  appear  to  be  closed  or  bounded.  Nevertheless,  this  problem  is  easily 
overcome  by  the  fact  that  once  we  "measure"  the  world  (i.e.,  collect  a  set  of 
data  points  experimentally)  we  reduce  the  real-world  phenomenon  to  a 
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dosed,  bounded  set  (since  we  can  only  collect  a  finite  number  of  pieces  of 
information).   Since  all  of  our  observations  of  nature  are  reduced  to  a  finite 
set,  the  closed  and  bounded  assumption  fits  nicely  with  our  view  of  the  real 
world. 

There  is,  however,  another  problem  with  the  experimental  method;  a 
problem  which  is  not  so  easily  dismissed.  Recall  that  the  definition  of  fractal 
dimension  invokes  the  limit  as  e  approaches  zero.   But  this  limiting  process 
cannot  be  reproduced  experimentally.  For  the  fractals  created  through 
infinite  iterations  of  function  systems,  the  attractors  are  known  down  to  the 
most  infinite  geometric  detail.  Moreover,  they  are  completely  self-similar  at 
every  stage  of  iteration.  However,  we  do  not  have  this  luxury  in  the  real 
world.  We  only  obtain  values  of  e  and   N(A,  e)  for  a  finite  set  of  nonzero  e 
values,  with  no  guarantee  that  the  observed  behavior  will  continue  as  e 
approaches  0.  For  example,  consider  a  piece  of  coral  which  appears  to  be 
similar  in  construction  to  a  fractal  tree.  On  the  macroscopic  level  we  find 
numerous  self-similar  aspects  of  the  geometric  shape  of  the  coral.  However, 
as  smaller  and  smaller  values  of  e  are  taken,  we  find  that  on  the  cellular 
level  there  is  no  similarity  between  the  shapes  of  the  cells  and  the  shape  of 
the  coral.   Further  analysis  reveals  that  the  cell  walls  themselves  even  have  a 
different  fractal  dimension.  Moreover,  this  problem  could  perpetuate 
endlessly  as  we  continue  to  decrease  e. 

As  another  example  of  this  problem,  consider  the  dimension  of  a  ball  of 
string.  From  a  distance,  the  ball  appears  to  occupy  three  dimensions  in  space, 
so  we  might  conclude  its  dimension  is  3.  However,  as  we  move  closer  to  the 
ball,  we  see  it  not  as  a  solid  object  but  as  a  long  string  wound  around  itself,  so 
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we  might  conclude  it  has  only  one  dimension.  As  we  get  closer  still,  we  see 
the  thickness  of  the  string,  and  again  we  conclude  it  is  three-dimensional,  but 
further  inspection  shows  many  one-dimensional  fibers  wound  together  to 
make  the  string.   These  examples  demonstrate  the  inherent  problem  with 
calculating  fractal  dimension  over  different  ranges  of  e,  and  we  have  not  yet 
taken  e  anywhere  near  0. 

Since  we  are  plotting  e  and  N(A,  e)  only  over  a  finite  range  of  e  values, 
any  conclusions  drawn  about  fractal  dimension  must  rely  on  the  gross 
assumption  that  the  object  is  infinitely  self-similar.   Consequently,  fractal 
dimension  has  been  more  often  used  to  disprove  hypotheses  than  to 
conclusively  prove  them.  Nevertheless,  if  two  objects  have  similar  fractal 
dimensions  over  the  same  range  of  e  values,  and  if  one  object  is  well 
understood,  one  may  use  this  information  to  help  explain  and  analyze  the 
other  object. 

As  an  example,  consider  the  dissipation  of  different  pollutants  in  a 
laboratory  "atmosphere."  In  an  isolated  setting,  the  behavior  can  be  studied 
in  great  detail  and  the  fractal  dimension  of  the  pollution  clouds  can  be 
determined  over  a  wide  range  of  e  values.  Observing  the  fractal  dimension  of 
real-world  pollution  clouds  over  the  same  values  of    e  might  be  a  way  to 
determine  which  pollutant  is  controlling  the  clouds'  behavior,  and  which 
pollutants  are  escaping  into  the  atmosphere.   Similar  techniques  have  been 
used  to  liken  the  growth  of  certain  plants  to  the  formation  of  crystals  in  a 
laboratory,  and  to  extrapolate  the  behavior  of  the  known  system  to  that  of  an 
unknown  system. 
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Two  other  techniques  for  experimentally  determining  fractal  dimension 
are  the  area-perimeter  relationship  and  the  number-resolution  relationship. 
In  the  area-perimeter  relationship,  the  area  of  a  set  of  data  is  computed,  and 
its  perimeter  is  measured  several  times  on  a  finer  and  finer  scale.   In  the 
number-resolution  relationship,  the  number  of  pixels  on  a  screen  is  plotted 
against  the  resolution  of  the  screen  for  a  variety  of  resolutions.  In  both  of 
these  techniques,  the  results  (graphs  of  ordered  pairs  of  perimeters  and  scale, 
or  numbers  and  resolutions)  are  usually  compared  to  those  of  a  known  set  of 
data  points.  If  the  slopes  of  these  data  points  (or  the  logs  of  the  data  points) 
correspond  to  an  object  with  a  known  fractal  dimension,  then  the  fractal 
dimension  of  the  unknown  object  can  be  estimated.   These  two  techniques 
are  used  widely  in  practice,  particularly  with  photographic  images,  but  again 
are  most  commonly  used  to  determine  the  similarity  between  two  fractal-like 
objects. 

Because  the  field  of  fractal  geometry  is  so  new,  the  possibilities  of 
applications  have  only  begun  to  be  discovered;  however,  fractal  dimension  is 
one  aspect  that  has  proven  quite  useful  in  many  sciences. 

J.     THE  KOCH  SNOWFLAKE 

An  interesting  phenomenon  that  has  recently  been  explained  using  these 
ideas  is  the  discrepancy  in  the  recorded  lengths  of  the  borders  of  many 
countries.  In  Europe,  the  reported  lengths  of  the  borders  between  many 
countries  varies  greatly,  depending  on  who  gives  the  report.  For  example,  the 
French  records  and  the  German  records  of  the  length  of  the  French-German 
border  may  differ  by  several  kilometers.  The  reason  for  this  difference  is 


73 


simple.  A  border  measured  with  a  ruler  from  a  1:50,000  scale  satellite 
photograph  records  much  less  of  the  border's  detail  than  a  surveying  team  on 
the  ground  records,  which  again  would  record  much  less  detail  than  a 
government  employee  on  the  ground  with  a  ruler.   The  smaller  the  device 
used  to  measure  the  border,  the  more  detail  will  be  recorded,  and  the  longer 
the  border  will  appear  to  be.  This  idea  is  readily  illustrated  through  the 
analysis  of  a  fractal  object  called  the  Koch  snowflake. 

To  construct  a  Koch  snowflake,  begin  with  an  equilateral  triangle,  and 
trisect  each  of  its  sides.  At  each  point  of  bisection,  construct  another 
equilateral  triangle  extending  outwards  from  the  original  triangle,  and  trisect 
each  of  their  external  sides  (see  Figure  3.28).  Continuing  this  geometric 
construction  indefinitely  yields  the  Koch  snowflake.   While  the  area  inside 
the  Koch  snowflake  is  clearly  finite  (since  the  entire  snowflake  lies  inside  any 
circle  containing  the  original  triangle)  the  length  of  its  border  is,  surprisingly, 
infinite.  To  see  this,  note  that  at  each  step  the  length  of  the  border  is 
increased  by  a  factor  of  4/3.  Thus,  starting  with  a  length  of  L  yields  a  final 
length  of  limn_>0oL(4/3)n,  which  is  an  infinite  length.  If  we  were  to  measure 
the  Koch  snowflake  with  rulers  of  different  scales,  remarkably  different 
results  for  the  length  of  the  border  will  be  obtained.  This  is  the  same 
phenomenon  being  experienced  measuring  border  lengths  in  Europe. 

As  another  real-world  application  of  this  concept,  many  small  "resort" 
lakes  entice  tourists  by  advertising  "a  thousand  miles  of  shoreline,"  or 
promises  to  that  effect.  Never  in  such  an  advertisement  will  one  find  a 
mention  of  the  scale  used  to  arrive  at  that  number,  as  this  is  simply  an 
application  of  the  above  concepts. 
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Second  Step 


Third  Step 
Figure  3.28  Constructing  the  Koch  snowflake. 

K.    APPLICATIONS  OF  FRACTAL  GEOMETRY 

Aside  from  the  applications  already  mentioned  (namely,  an  aid  to 
computer  graphics  and  the  use  of  fractal  dimension  to  classify  and  compare 
similar  objects)  fractal  geometry  has  so  far  found  few  other  applications  in  the 
physical  sciences.  As  already  mentioned,  the  infinitely-detailed  property  of 
fractals  does  not  hold  in  the  real  world  when  objects  are  analyzed  at  the 
molecular  and  smaller  levels.   Even  Benoit  Mandelbrot,  the  German-born 
American  scientist  who  gave  fractals  their  name,  admits  that  true  fractals  do 
not  exist  in  nature.  However,  he  is  also  quick  to  point  out  that  there  are  no 
truly  straight  lines  or  perfect  circles  either.  While  traditional  geometries  have 
been  inadequate  in  sufficiently  describing  all  of  nature,  perhaps  fractal 
geometry,  at  least  on  the  macroscopic  level,  may  provide  yet  another 
approach.  When  we  consider  the  infinitude  of  natural  objects  that  have  a 
fractal-like  structure,  we  see  that  nature  is  indeed  more  fractal  than  it  is 
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Euclidean.  From  the  formation  of  crystals  and  snowflakes  to  the  way  coral, 
certain  roots,  plants,  and  trees  grow;  from  the  shapes  of  lightning  bolts  and 
some  electrical  discharges  to  cloud  formations,  weather  patterns,  and  galactic 
patterns;  from  the  structure  of  our  lungs  and  surfaces  of  our  brains,  to  the 
patterns  of  our  veins;  much  of  the  physical  world  seems  to  have  a  fractal 
structure  to  it.  A  more  complete  analysis  of  some  of  the  minutely  detailed, 
self-similar  objects  in  nature  may  someday  be  realized  through  fractal 
geometry. 

One  area  of  mathematics,  however,  where  fractal  geometry  has 
unquestionably  found  a  permanent  role  and  distinguished  itself  is  in  the 
study  of  chaotic  dynamical  systems.  Dynamical  (changing)  systems  are 
prevalent  in  the  real  world  and  we  are  finding  ever  greater  numbers  of  them 
to  be  chaotic.  In  both  mathematical  models  and  the  physical  world,  strange 
sets  of  data  points  are  being  observed  that  seem  to  contain  self-similarity  and 
infinite  detail.  Many  of  these  sets  of  points  are,  in  fact,  fractals  and  with  the 
aid  of  fractal  geometry  we  can  understand  them  better.  Hence,  when  we 
realize  how  much  of  our  universe  is  modeled  by  chaotic  dynamical  systems, 
we  will  better  appreciate  the  usefulness  of  fractal  geometry  in  studying  these 
systems.  We  take  up  the  study  of  chaos  in  the  next  chapter. 
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IV.      CHAOS 

A.    INTRODUCTION 

We  have  introduced  and  studied  iterated  function  systems,  from  our 
simplest  example  where  fn+1(x)  =  (PKx))2,  to  systems  that  produced  very 
intricate  and  interesting  geometric  shapes.  The  similarity  between  these 
systems  is  that  we  could  predict  the  exact  behavior  of  the  systems  after  any 
number  of  iterations.  Unfortunately/  in  nature  this  predictability  seems  to  be 
more  the  exception  than  the  rule,  which  has  caused  scientists  frustration  for 
centuries  whenever  their  models  of  seemingly  simple  phenomena  produced 
erratic  behavior.  This  observed  phenomena  in  dynamical  systems  is  called 
chaos,  a  precise  definition  of  which  we  provide  below. 

We  introduce  our  study  of  chaotic  dynamical  systems  by  first  looking  at 
discrete  systems,  such  as  the  iterated  function  systems  we  have  been 
discussing.  There  are  two  reasons  for  this  approach:  the  first  is  that  discrete 
systems  provide  a  conceptually  simpler  setting  for  understanding  the  theory, 
and  the  second  is  because  many  data  are  collected  from  the  real  world  in 
discrete  increments.  For  example,  consider  a  biologist  studying  the  dynamics 
of  a  population.  Not  only  will  the  population  change  in  discrete  increments, 
but  the  biologist  can  measure  the  population  only  at  discrete  intervals.  In 
fact,  even  in  continuous  real-world  behaviors  such  as  beam  vibrations,  we  are 
able  to  record  measurements  only  at  discrete  intervals.  Since  much  of  science 
involves  studying  a  continuous  real  world  based  on  discrete  observations,  it 
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makes  sense  to  begin  the  study  of  chaotic  dynamical  systems  with  the  discrete 
case. 

A  further  reason  to  study  discrete  iterated  function  systems  is  that  they 
capture  the  essence  of  an  important  real-world  phenomenon  known  as 
feedback.  Feedback  occurs  when  the  present  state  of  a  system  affects  its  future 
state,  which  in  turn  means  there  could  be  at  least  a  small  time  delay  in  the 
change  of  the  state  of  a  system.  For  instance,  feedback  in  speaker-microphone 
systems  causes  jumps  at  discrete  time  intervals  because  of  the  time  it  takes  the 
sound  to  travel  through  the  system  (although  we  preceive  these  jumps  as 
continuous  when  we  hear  them).   Biological  systems  which  are  periodic  and 
demonstrate  feedback  include  predator-prey  models,  the  motion  of  slime 
molds,  light  emission  by  groups  of  fireflies,  glycolysis  and  photosynthesis,  and 
even  the  conditions  in  the  brain  leading  to  epileptic  seizures.  Feedback  is  also 
exhibited  in  electrical  circuits,  as  the  defense  contractor  TRW  discovered 
when  chaotic  feedback  shut  down  its  European  computer  network.  Finally, 
gaps  in  the  asteroid  belts  that  would  be  in  phase  with  planetary  orbits  indicate 
the  effect  of  feedback  in  our  solar  system.  These  are  just  a  few  examples 
where  relatively  simple  mathematical  models  can  capture  important  real- 
world  phenomena. 

The  following  five  sections  are  based  on  Devaney  (1989),  and  provide  a 
background  and  a  setting  for  studying  chaotic  dynamical  systems.  The 
remainder  of  the  chaper  comes  from  a  variety  of  sources,  but  in  the  cases 
where  the  references  differ  on  theorems  or  definitions,  the  ones  used  by 
Devaney  are  followed.  Additionally,  many  of  the  examples  and  applications, 
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where  not  specifically  referenced,  come  from  Gleick  (1988),  Moon  (1987),  Levy 
(1991),  and  Briggs  (1989). 

We  begin  by  developing  the  mathematics  of  discrete  dynamical  systems. 
Our  goal  is  to  understand  thoroughly  the  idea  of  orbits  of  points  under 
repeated  iteration. 

B.    GRAPHICAL  ANALYSIS  OF  FIXED  POINTS  OF  MAPS 

We  have  already  shown  how  graphical  analysis  of  a  function  can  reveal 
information  about  its  fixed  points,  periodic  points,  and  whether  those  points 
are  attracting  or  repelling.  We  will  now  introduce  the  concept  of  graphical 
analysis  of  the  nth  iterate  of  a  function,  denoted  by  f^x).  Consider  first  a  very 
simple  function  mapping  the  unit  interval  to  itself   f:  [0,  l]-»[0, 1],  defined  by 
f(x)  =  2x  (modi).  The  notation  2x  (modi)  means:  multiply  the  input  value 
by  two,  and  consider  only  the  fractional  part,  i.e.,  the  remainder  upon 
division  by  one.  For  example,  f(1.57)  =  3.14  (modi)  *  .14.  This  mapping  gives 
f(x)  =  2x  if  0  <,  x  £  1/2  and   f(x)  m  2x  - 1  if     1/2  <  x  £  1.  This  function  is  called 
the  Baker  map,  and  its  graph  is  shown  in  Figure  4.1.  The  Baker  map  is  clearly 
an  onto  function,  but  it  is  not  one-to-one  since  every  y  value  in  [0, 1]  is  the 
image  of  exactly  two  values  of  x.  Additionally,  the  fixed  points  of  the  Baker 
map  are  0  and  1,  both  of  which  are  repelling. 

Consider  the  second  iterate  of  the  Baker  map,  f2(x)  =  4x  (modi).  Its 
behavior  can  be  seen  more  clearly  by  the  graphical  analysis  shown  in  Figure 
4.2.  Since  the  interval  [0,1/2]  is  mapped  to  [0,1]  under  the  first  iterate,  it 
makes  sense  that  under  the  second  iterate  the  interval  [0, 1/2]  will  contain  a 
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►    X 


u  1/2  1 

Figure  4.1  The  Baker  map. 

smaller  copy  of  the  entire  first  map.  Additionally,  the  second  iterate  shows 
the  appearance  of  two  new  fixed  points,  1/3  and  2/3.    This  fact  demonstrates 
the  existence  of  a  two-cycle  in  the  original  map;  in  fact,  it  is  easy  to  verify  that 
f(l/3)  =  2/3  and   f(2/3)  - 1/3.  Notice  also  that  the  fixed  points  of  f(x),  0  and 
1,  remain  fixed  for  f2(x).  Observing  the  fixed  points  of  the  nth  iterate  of  a 
map  reveals  both  the  periodic  and  fixed  points  of  the  original  map.  The  third 
iteration  is  f3(x)  =  8x  (modi)  and  its  fixed  points,  0, 1/7, 2/7, . . . ,  6/7, 1    are 
also  shown  in  Figure  4.2.  In  the  case  of  the  Baker  map,  we  can  easily  obtain  a 
closed-form  expression  for  its  nth  iterate:    f^x)  =  2nx  (modi).  This  formula 
makes  the  Baker  map  very  easy  to  graph  and  analyze,  and  we  will  return  to  it 
when  "chaos"  is  defined.  It  may  be  very  difficult,  however,  to  obtain  a  simple 
expression  for  the  nth  iterate  of  a  map,  and  it  is  frequently  even  more 
difficult  to  graph  it. 

A  slightly  more  complicated  map  to  analyze  graphically  is  the  map 
f(x)  =  4x(l  -  x),  which  is  a  parabola  that  maps  the  unit  interval  to  itself.  Its 
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graph  and  second  iterate  on  the  unit  interval  are  shown  in  Figure  4.4.   Again, 
the  appearance  of  fixed  points  under  f2(x)  shows  a  two-cycle  of  f(x),  and  a 
simple  check  shows  this  two-cycle  occurs  at  the  points  x  =  (5  ±  V5)/8. 


f     (x)  =  4x(mod1)  f     (x)  =  8x(mod1) 

Figure  4.2  The  second  and  third  iterates  of  the  Baker  map. 


1/2  1*X 

f(x)  =  4x(1  -  x) 
Figure  4.3  The  first  and  second  iterates  of  f(x)  =  4x(l  -  x). 


►  x 
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While  this  technique  is  not  always  easy  to  use  (due  to  the  difficulty  of 
graphing  many  functions  accurately),  it  can  reveal  considerable  information 
about  the  orbits  of  points  in  a  dynamical  system. 

C    MAPS  OF  THE  CIRCLE 

Consider  the  map  f:  S1-^1  with  f(9)  =  26  of  the  unit  circle  to  itself.  Here 
0  is  the  angle  measured  in  radians  (positive  counterclockwise)  between  the 
positive  x-axis  and  the  line  joining  the  origin  to  any  point  on  S1.  Every  point 
on  the  circle  is  represented  by  the  form  050  +  2krc  for  any  integer  k.  The 
map  f,  which  doubles  the  angle  of  any  point  on  the  unit  circle,  has  certain 
properties  which  help  illustrate  and  suggest  a  precise  definition  of  chaotic 
systems  (see  Figure  4.4). 


Figure  4.4  The  map  f(0)  =  20  of  the  circle  S1. 

First  notice  that  0  =  0  is  a  fixed  point  of  f,  since  f(0)  =  0.  Additionally,  for 
any  point  of  the  form  0  =  2k7t/2n,  ^(0)  =  2krc,  so  that  these  points  are 
eventually  fixed.  The  points  of  the  form  0  =  2krc/2n  form  a  dense  subset  of 
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S1,  since  between  any  two  distinct  points  on  S1,  we  can  find  a  point  of  the 
form  2krc/2n.  To  see  this,  consider  two  arbitrary  distinct  points  9i  and  82  on 
S1,  where  61  =  2na\  and  02  =  2na.2  with  0  <  ai  <  0C2  <  1.  Since  02  -  ai  >  0, 
there  exists  an  integer  m  such  that  m((X2  -  ai)  >  1  by  the  Archimedian 
property.  Let  2n  >  m.  Then  2n(a2  -  ai)  >1  implies  there  exists  a  positive 
integer  k  such  that  2noi2  >  k>  2nai,  or  62  >  2krc/2n  >  61.  Thus,  points  of  the 
form  2krc/2n  are  dense  on  S1. 

The  periodic  points  of  f  are  slightly  more  difficult  to  find.  Since 
1^(0)  =  2n0,  the  point  0  is  periodic  of  period  n  if  and  only  if   2n0  =  0  +  2krc 
for  some  integer  k.  Solving  for  0  yields  0(2n  -  1)  =  2krc,  so  for   0  <,  k  £  2n  -  1, 
the  points  0  =  2krc/(2n  -  1)  are  the  periodic  points  of  period  n.  Using  complex 
analysis,  it  can  be  shown  that  these  points  are  the  (2n  -  1)  roots  of  unity,  and 
that  they  are  dense  on  the  circle  S1.  Geometrically,  this  means  that  the  points 
of  period  n  are  those  that  divide  the  unit  circle  into    2n  -  1  equal  segments. 
This  interpretation  also  shows  that  they  are  dense  on  S1.  This  is  shown 
graphically  in  Figure  4.5  for  the  points  of  period  two.  The  point  0o  trivially 
has  period  two  since  PKOo)  =  0o  for  all  n.  The  points  0i  and  02  have  period 
two  since  f(0i)  =  02  and  f(02>  =  H 

Another  interesting  property  of  f  is  that  for  any  two  points  on  S1,  no 
matter  how  close  together  they  start,  iteration  under  f  will  separate  their 
iterates  by  an  arbitrary  amount.  The  formal  definition  of  this  characteristic  is 
as  follows:  a  mapping  f:  S-*S   has  sensitive  dependence  on  initial  conditions 
if  there  exists  5  >  0  such  that  for  any  point  x  e  S  and  any  neighborhood  N 
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e  = 


92  =    4tc/3 
Figure  4.5  Periodic  points  on  S1  of  period  two. 


of  x,  there  exists  a  point  y  €  N  and  an  integer  n  >  0  such  that 
d(fn(x),  PKy))  >  8.  Notice  that  under  this  definition,  it  is  not  necessary  for 
every  point  in  the  neighborhood  N  to  move  away  from  x  under  repeated 
iteration.   However,  in  any  neighborhood  containing  x  there  must  exist  at 
least  one  point  that  does  move  away  from  it.  Since  our  map  of  the  circle 
f(6)  =  20  doubles  the  aro  length  between  any  two  points,  there  is  always  an 
iteration  number  n  that  moves  two  points  at  least  a  given  distance  8  apart; 
hence  the  map  f(0)  =  20   has  sensitive  dependence  on  initial  conditions. 

A  third  interesting  property  of  the  map    f(0)  =  20  is  that  for  any  pair  of 
open  arcs  U,  V  c  S1,  there  exists  a  number  n  >  0  such  that  the  intersection  of 
f^flJ)  and  V  is  nonempty.  The  formal  definition  of  this  property  is  as 
follows:  a  mapping  f:  S-»S  is  called  topologically  transitive  if  for  any  pair  of 
open  sets  U,  V  c  S  there  exists  an  integer  n  >  0  such  that    f^CU)  nV^0. 
Since  any  arc  U    in  S1  is  expanded  under  iteration  of  f(0)  ■  20  until  it 
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covers  all  of  S1,  every  arc  eventually  intersects  every  other  arc  in  S1.  Thus 
the  mapping  f(9)  =  29  is  topologically  transitive. 

The  three  properties  of  f:  S1-^1  where  f(9)  =  29  define  a  dynamical 
system  to  be  chaotic.   The  formal  definition  follows. 

D.    CHAOTIC  DYNAMICAL  SYSTEMS 

A  mapping  f:  S-*S  is  called  chaotic  if  the  following  three  conditions  are 
met: 

1.  f  has  sensitive  dependence  on  initial  conditions; 

2.  f  is  topologically  transitive;  and 

3.  Periodic  points  are  dense  in  S. 

While  precise  definitions  of  chaos  vary  slightly  among  different 
disciplines,  most  of  them  are  at  least  similar  to  those  we  gave  here.  These 
three  conditions  are  frequently  referred  to  as  unpredictability, 
indecomposability,  and  an  element  of  regularity,  respectively.  Because  of  the 
sensitive  dependence  on  initial  conditions,  it  is  almost  impossible  to  predict 
the  orbit  of  an  arbitrary  point  with  any  degree  of  accuracy.  This  is  the 
unpredictability  aspect  of  chaos,  and  it  makes  numerical  computations  with 
chaotic  dynamical  systems  virtually  impossible  (since  the  slightest  rounding 
error  at  any  step  will  almost  certainly  create  a  point  with  an  entirely  different 
orbit).  The  indecomposability  property  comes  from  the  fact  that  there  is  no 
way  to  decompose  the  set  S  into  two  disjoint  subsets  that  do  not  eventually 
interact  (since  any  open  subset  of  S  eventually  intersects  every  other  open 
subset  of  S  under  iteration  of  f).  Finally,  the  denseness  of  the  periodic  points 
provides  an  element  of  regularity  amidst  all  this  chaotic  behavior. 
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Let  us  now  demonstrate  that  the  Baker  map  is  chaotic  directly  through  a 
graphical  analysis.  To  see  that  the  periodic  points  are  dense,  note  that  the  nth 
iterate  f^Kx)  =  2nx  (modi)  has  n  evenly  spaced  fixed  points  in  [0, 1]  for  n 
arbitrarily  large  as  indicated  in  Figure  4.2.  Figure  4.6  also  shows  that  f  is 
topologically  transitive  since  any  open  interval   (p,  q)   contains  an  interval 
that  is  mapped  to  [0,1].  Finally,  Figure  4.7  shows  the  sensitive  dependence 
on  initial  conditions,  because  any  two  points  p  and  q  are  moved  an  arbitrary 
distance  apart  under  some  iteration  of  f.  While  this  graphical  technique  is 
not  a  proof,  it  can  be  supported  by  a  rigorous  mathematical  argument.  It 
alsodemonstrates  the  usefulness  of  analyzing  some  maps  graphically  to 
determine  if  they  are  chaotic. 


f(x) 


■f 


+ 


p  q  x    (not  to  scale) 

Figure  4.6  Topological  transitivity  of  the  Baker  map. 
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f(q) 


Af(x) 


f(p) 


p         q  x    (not  to  scale) 

Figure  4.7  Sensitive  dependence  on  initial  conditions  of  the  Baker  map. 

Chaos  occurs  not  only  in  mathematical  equations,  but  in  fact  seems  to  be 
most  prevalent  in  the  real  world.   While  accurately  modeling  real-world 
phenomena  with  mathematical  equations  is  sometimes  very  difficult, 
particularly  when  the  physical  system  is  chaotic,  it  happens  that  even  a 
simple  model  crudely  representing  a  physical  system  may  turn  out  to  be 
chaotic.  This  happenstance  in  turn  tells  us  that  the  more  complex  physical 
system  is  probably  also  chaotic.  We  will  see  this  happen  when  we  discuss  the 
Lorenz  equations  as  a  crude  model  for  weather  prediction.  Another  example 
of  this  chaotic  behavior  occurs  with  a  simple  differential  equations  model 
used  to  describe  beam  intersections  in  particle  accelerators.  The  model 
explains  why  it  is  so  difficult  to  predict  the  action  of  intersecting  beams  in 
particle  accelerators  (at  least  using  existing  models). 

There  are  several  examples  of  physical  systems  where  chaos  has  been 
observed,*  studied,  and  in  some  cases,  somewhat  understood.  These  systems 
include  turbulence  in  fluids,  thermal  convection  in  gasses,  panel  flutter  on 
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supersonic  aircraft,  certain  chemical  reactions  (specifically  the  Belousof- 
Zhabotinsky  reaction),  abnormal  cardiac  rhythms,  nonlinear  electrical 
circuits,  biological  population  dynamics,  vibrations  of  buckled  elastic  systems 
(such  as  beams),  geomagnetic  field  reversals,  and  even  planetary  motion. 
Many  of  these  systems  are  discussed  in  Moon  (1987),  Holden  (1986),  and 
Rasband  (1990). 

Unfortunately,  it  is  very  difficult  to  identify  chaotic  systems 
experimentally.  One  reason  is  that  numerical  roundoff  in  collected  data 
could  lead  to  an  erroneous  assumption  that  two  data  points  are  the  same 
(indicating  a  cycle),  whereas  a  slight  difference  in  their  actual  values  may 
cause  their  orbits  to  diverge  rapidly  (because  of  sensitive  dependence  on 
initial  conditions).   An  even  more  difficult  problem  is  with  distinguishing 
between  a  truly  chaotic  orbit  and  a  cycle  with  a  very  long  period.  To 
emphasize  the  significance  of  this  problem  we  note  that  the  Department  of 
Mechanics  at  Cornell  University  requires  4,000  non-cyclic  real-world  data 
points  before  scientists  there  declare  a  system  as  being  chaotic.  For  biologists 
and  economists,  who  model  annual  trends,  this  requirement  may  seem 
unreasonable.  However,  to  a  helicopter  pilot  whose  life  depends  on  chaotic 
vibrations  not  occurring  on  the  main  rotor  blade,  this  requirement  may  not 
be  stringent  enough. 

While  the  map  f:  S1-^1  with  f(9)  =  29  provides  a  simple  example  of  a 
chaotic  map,  it  by  no  means  represents  the  extent  of  chaotic  dynamical 
systems.  We  turn  now  to  several  more  interesting  chaotic  maps. 
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E    TOPOLOGICAL  CONJUGACY 

One  way  to  determine  if  a  map  is  chaotic  is  to  check  directly  the  three 
conditions  of  chaos.  However,  using  the  concept  of  topological  conjugacy,  we 
can  apply  the  dynamics  of  a  familiar  map  to  those  of  another  function.  Two 
maps,  f :  X-»X  and  g:  S-»S  are  said  to  be  topologically  conjugate  if  there  exists 
a  homeomorphism  h:  S-»X  such  that  f(h)  =  h(g).  If  f  and  g  are 
topologically  conjugate,  then  g  =  h~*fh,  so  gn  =  (h-1fh)n  =  (h^fhXh^fh) . . . 
(h^fh)  =  (h-1f)(hh-1)f(hh-1) . . .  fOuYtyfh).  Since  h  is  a  homeomorphism, 
hh_1  is  the  identity  function,  hence  gn  =  h'^h,  and  hgn  =  f^h.  It  follows  that 
f  and  g  share  the  same  dynamical  properties;  in  particular  if  g  is  chaotic  the 
same  is  true  of  f.  This  idea  is  shown  in  Figure  4.8  with.  If  h:  S-»X  is  not  a 
homeomorphism  (e.g.,  if  it  is  two-to-one)  then  f  and  g  are  said  to  be 
topologically  semi-conjugate.  Nevertheless,  if  g  is  chaotic,  f  is  still  also 
chaotic.  Since  h  is  two-to-one,  the  dynamics  of  f  are  even  more  complicated 
because  h  introduces  even  more  periodic  points  or  cycles.  It  suffices  that 
chaotic  behavior  is  preserved  through  topological  semi-conjugacy,  which  is 
the  extent  to  which  we  use  this  concept. 

S      —  g   ->     S 

in  in 

X      —  f    ->    X 

Figure  4.8  Diagram  showing  topological  conjugacy  between  g  and  f . 

What  follows  is  an  example  of  the  usefulness  of  topological  conjugacy. 
Consider  the  map  f:  R-»R,  where  f(x)  =  2x2  -  1.  We  would  like  to  know  if 
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this  map  is  chaotic.  We  could  try  to  show  sensitive  dependence  on  initial 
conditions,  topological  transitivity,  and  search  for  a  dense  orbit.  Instead  we 
consider  the  mapping  h(9)  =  cos9.  Since,  cos26  =  2cos28  - 1  =  f(h(0)),  the  map 
f  is  topologically  semi-conjugate  (since  the  cosine  function  is  not  one-to-one) 
with  the  map  g(9)  =  20  when  g:  S1-»S1'  as  shown  in  Figure  4.9.  Since  we 
know  that  g  is  a  chaotic  map  from  our  previous  discussions,  and  that  f 
shares  the  same  dynamical  properties  as  g,  the  map  f(x)  =  2x2  -  1  is  also 
chaotic. 


9      —  g    ->     29 
ih  ih 

cos0  —  f    — »    cos29 


Figure  4.9  Topological  conjugacy  between  g(9)  =  29  and  f(x)  =  2x2  -  1, 
under  h(9)  =  cos9. 


We  frequently  use  the  concept  of  topological  conjugacy  to  analyze  maps 
suspected  to  be  chaotic.  The  simplicity  of  the  previous  example  shows  why 
topological  conjugacy  is  often  used. 

F.    CHAOTIC  DYNAMICS  ON  CODE  SPACE 

Consider  now  the  shift  map  a:  £2-^X2  on  code  space    Z2  =  {S1S2S3. . . : 

Si  €  {0, 1}}  where  o(siS2S3- . .)  ■  S2S3S4 We  show  directly  from  the 

definition  that  this  mapping  is  chaotic. 

To  see  that  a:  Z2~>22   has  sensitive  dependence  on  initial  conditions, 
consider  any  point  x  =  X1X2X3. . .  in  £2-  For  a  point  y  =  yiy2y3- . .  to  be  in  a 
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neighborhood  of  x,  it  must  agree  with  x  for  a  certain  number  of  digits 
depending  on  the  size  of  the  neighborhood.  Suppose  the  agreement  is  with 
the  first  n  digits;  xi  =  yi,  X2  =  y%  •  •  •  /  xn  -  Yn-  If  we  choose  yk  *  x^   for  all 
k>  n,  then  y  still  lies  in  the  neighborhood  of  x,  but   I  on(x)  -  on(y)  I  >  8  for  all 
8  <  1.  Hence,  the  point  y,  which  is  arbitrarily  close  to  x,  moves  arbitrarily 
away  from  it  under  iteration  of    a.  This  establishes  sensitive  dependence  on 
initial  conditions. 

To  show  topological  transitivity,  consider  two  open  sets  U,  V  €  Yji-  Since 
U  is  open,  its  points  U1U2U3. .  .un. . .  in  £2  agree  for  the  first  n  digits,  but 
can  differ  in  any  way  beyond  the  nth  digit  for  some  positive  integer  n. 
Likewise,  points  in  £2  which  agree  in  the  first  m  digits  for  some  m  all  lie  in 
V.  It  follows  that  the  point  U1U2U3. .  .unviV2V3. .  .vm. . .  lies  in  U,  while  the 
nth  iterate  of  U  under  a  contains  the  point  V1V2V3. .  .vm. . .  in  V.  Thus  the 
intersection  of  an(U)  and  V  is  nonempty  for  any  open  sets  U  and  V  of  £2- 
Therefore,  a:  Y/i-^Li  is  topologically  transitive. 

To  see  that  the  periodic  points  of  a  are  dense  in  X2/ w^  show  that  for  any 
point   se  I2  there  is  a  periodic  point  arbitrarily  close  to  it.  Thus,  let 

s  =  S1S2S3 We  want  to  find  a  periodic  point   x  that  corresponds  to  s  for 

up  to  n  digits  (selecting  n  arbitrarily  large  makes  the  periodic  point 
arbitrarily  close  to  s).  Choose  the  point  x  =  S1S2. .  .SnSiS2- .  .Sn- . . ,  which  is  a 
periodic  point  of  period  n.  The  point  x  is  arbitrarily  close  to  S.  Hence,  the 
periodic  points  of  o  are  dense  in  £2- 

In  conclusion,  the  shift  map  on  code  space  a:  X2->£2  is  chaotic.  Like  the 
map  of  the  unit  circle  f(0)  =  28,  the  shift  map  is  quite  useful  in  showing  that 
other  maps  are  chaotic  through  the  use  of  topological  conjugacy.  When  we 
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can  establish  topological  conjugacy  between  the  maps  f  and  o,  we  will  refer 
to  analyzing  the  behavior  of  points  under  f  as  studying  the  symbolic 
dynamics  of  the  map  f . 

Using  the  concept  of  topological  conjugacy,  it  is  now  an  easy  exercise  to 
show  that  the  Baker  map  f(x)  =  2x  (modi)  is  chaotic.  This  behavior  can  be 
shown  through  either  the  shift  map  on  code  space  or  the  angle  doubling  map 
on  the  unit  circle,  which  confirms  our  earlier  geometric  demonstration  of 
this  fact. 

a    NEWTON'S  METHOD  FOR  X2  =  -1 

To  see  how  chaos  sometimes  arises  from  mathematical  systems,  consider 
the  following  example  from  numerical  analysis  due  to  Strang  (1991).  Recall 
that  in  finding  the  zeros  of  a  real  valued  function  f(x)  by  Newton's  method, 
an  initial  guess  together  with  the  iterative  process  xn+i  -  xn  -  f(xn)/f  (xn)  is 
used.  Applying  this  method  to  find  the  roots  of  f(x)  =  x2  +  1,  the  iteration 
equation  reduces  to  the  iterated  function  system  xn+i  =  xn  -  (xn2  +  l)/(2xn)  = 
(l/2)(xn  -  1/xn).  Newton's  method  converges  for  most  polynomials  under 
certain  conditions  (e.g.,  when  f  (x)  *  0  and  no  inflection  point  occurs  between 
the  root  and  the  initial  guess).  In  the  case  in  question,  the  only  roots  are  ±i 
so  the  system  cannot  converge  on  the  real  line.  (But  note  that  f (x)  is  always 
strictly  positive  and  has  no  real  roots.) 

The  behavior  of  the  iterated  function  system  changes  dramatically 
depending  on  the  initial  guess.  For  some  values  of  xo,  this  system  diverges  to 
infinity;  for  example,  if  xo  =  1,  then  xi  =  0,  and  x2-»«»  (since  we  cannot 
divide  by  zero).  On  the  other  hand,  if  we  choose  xq  to  be  very  large,  then  xi 
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is  approximately  half  that  value,  and  the  orbit  moves  towards  zero  until  the 

l/xn  term  makes  the  next  iteration  large  again.  Finally,  if  we  choose 

xo  =  1/V3,  then  xi  =  -1/V3  and  X2  =  1/^3,  yielding  a  two  cycle. 

Because  we  are  searching  for  a  root  value  which  is  nonexistent  on  the 

real  line,  Newton's  method  exhibits  strange  behavior  for  x2  +  1  =  0.  We  now 

explore  whether  this  system  is  chaotic.  From  trigonometry,  cot20  = 

(l/2)(cot9  -  1/cote).  Using  the  map  h(0)  =  cotO  it  follows  that  the  map 

g(9)  =  29  is  topologically  conjugate  to  f(x)  =  (l/2)(x  -  1/x),  as  shown  in  Figure 

4.10. 

6      —  g   -»     29 
i  h  ih 

cot9  —  f   ->    cot29 


Figure  4.10  Topological  conjugacy  between  g(9)  =  29  and  f(x)  =  (l/2)(x  -  1/x) 
under  h(9)  =  cot9. 


Thus  f(x)  exhibits  the  same  dynamical  properties  as  g(9)  ■  29  on  S1;  hence  f 
is  chaotic.  We  can  actually  observe  this  behavior  graphically  by  considering 
the  graph  of  the  iterates  of  f .  Letting  xo  =  cot9,  we  obtain  xi  =  f  (xo)  =  cot29, 
X2  =  cot49,  and  in  general,  xn  =  cot(2n9).  Hence,  the  orbit  of  x  under  f  is  the 
sequence  cot9,  cot29,  cot49. . .  for  different  values  of  9  corresponding  to  our 
initial  guess.  The  graph  of  cot9  is  shown  in  Figure  4.11,  and  helps  to 
demonstrate  the  results  noted  earlier  for  different  initial  guesses. 

If  we  start  with  xo  =  1  we  get  9o  =  tc/4  (since  cot  (jc/4)  =  1).  Then  9i  =  rc/2, 
corresponding  to  xi  -  0  (since  cot  (rc/2)  =  0).  Next  92  =  Jc  and  the  orbit  for  x 
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f(G) 


►  e 


Figure  4.11  The  cotangent  map  and  a  chaotic  orbit. 


diverges  since  cot  n  is  undefined.  In  fact,  any  starting  value  of  xo 
corresponding  to  the  form  Go  =  k7c/2n,  for  positive  integers  k  and  n,  results 
in  an  orbit  that  eventually  diverges  to  ~.  Starting  with  xo  =  1/VJ,  we  find 
00  =  n/3  (since  cotOt/3)  =  1/V3).  Then  Gi  =  29o  =  2tc/3  gives  xi  =  -1/V5  (since 
cot(2;c/3)  =  -1/V3).  Next  62  =  29i  -  4n/3  =  Go  +  n.  Since  the  cotangent  function 
is  periodic  with  period  it,  X2  =  1/V3.  Therefore,  this  orbit  is  the  two-cycle  we 
observed  earlier.  In  fact,  if  Go  =  (p/q)rc,  where  p/q  is  not  of  the  form  k/2n, 
then  the  orbit  eventually  cycles.   This  observation  further  demonstrates  the 
denseness  property  of  the  periodic  points.  Finally,  if  60  is  an  irrational 
multiple  of  n,  the  orbit  will  never  diverge  or  cycle;  in  fact,  it  will  be  chaotic,  as 
we  already  knew. 

Considering  our  search  for  i  on  the  real  line,  the  points  that  we  hope 
would  eventually  converge  to  zero  are  yn  =  xn2  +  1.  If  we  further  analyze 
these  points,  we  find  yet  another  chaotic  map  which  shares  the  dynamics  of 
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the  one  we  just  studied.  For  yn+i  =  (xn+i)2  +  1  =(l/4)(xn  -  l/xn)2  + 1  = 
(l/4)(xn2  +  2  +  l/x„2).  Simplifying  algebraically,   yn+i  =  (l/4)(xn2  +  D2/xn2  = 
yn2/4(yn  -  1).  If  we  change  variables  and  let  z  =  1/y,  this  last  equation  reduces 
to  zn+i  =  4zn(l  -  zn).  The  latter  iterated  function  is  a  member  of  the  quadratic 
family  of  maps  we  study  in  great  detail  below,  and  we  know  it  is  chaotic. 

As  a  final  illustration  of  Newton's  method,  we  analyze  the  system 
f(z)  =  z4  -  1  in  the  complex  plane.  Since  this  equation  has  four  roots,  ±1  and 
±i,  depending  on  the  initial  value,  the  iterations  could  end  up  at  any  one  of 
the  four  roots.  The  results  encountered  with  computers  when  searching  for 
these  roots  has  not  always  been  predictable,  depending  on  the  initial  value. 
Now  each  of  these  roots  has  around  it  a  "basin  of  attraction"  inside  of  which 
all  initial  values  converge  to  that  root  upon  repeated  iteration.   However, 
along  the  diagonals  Re(z)  =  Im(z)  and  Re(z)  =  -Im(z),  there  are  regions  that 
give  rise  to  chaotic  behavior.  To  see  this,  color  the  complex  plane  based  on 
which  initial  values  converge  to  which  root:  that  is  color  all  the  points  that 
converge  to  1  blue,  to  -1  green,  to  i  yellow,  and  to  -i  red.  There  are  four 
large  basins  of  attraction  with  infinitely  detailed,  multicolored  borders 
between  them.  In  fact,  there  is  no  clearly  defined  boundary  between  any  two 
basins  of  attraction:  between  any  blue  and  red  region,  there  is  a  yellow  and  a 
green  one.  A  schematic  of  the  complex  plane  showing  the  basins  of  attraction 
and  the  chaotic  regions  is  shown  in  Figure  4.12.  A  detailed  color  image  of  this 
figure  can  be  found  in  Gleick  (1988,  p.  114). 
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Im(z) 


Chaotic 

Regions 
(Shaded) 


►  Re(z) 


Figure  4.12  Schematic  of  the  complex  plane  indicating  the  behavior 
under  Newton's  method  for  z4  +  1  and  the  (colored)  basins  of  attraction. 


H.   THE  QUADRATIC  FAMILY  OF  MAPS 

In  this  section  we  study  a  specific  case  of  maps  from  the  quadratic  family 
f(x)  =  ax2  +  bx  +  c.  We  discuss  a  model  from  population  biology  to  develop 
this  map,  which  comes  from  Briggs  (1989).  The  mathematics  in  this  section  is 
from  the  book  by  Devaney  (1989)  and  the  article  by  Devaney  (1989),  while 
some  of  the  more  detailed  results,  particularly  about  the  Feigenbaum 
constant,  are  from  Rasband  (1990). 

Consider  a  population  of  moths  which  live  for  one  season,  lay  their  eggs, 
and  then  die,  all  in  a  fixed  and  limited  environment.  If  we  do  not  consider 
the  interaction  of  the  moths  with  the  environment  or  other  species,  we  can 
model  the  dynamics  of  the  moth  population  by  assuming  that  the  size  of  the 
population  increases  after  each  life  cycle  in  constant  proportion  to  the  size  of 
the  population  during  the  previous  cycle.  If  we  let  the  constant  of 
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proportionality  equal  B  (for  birthrate),  then  for  each  population  cycle  we 
have  xn+i  =  Bxn/  where  x  represents  the  population  and  n  represents  the 
time  period.   Given  only  these  conditions  together  with  B  >  1,  this  model 
predicts  the  population  increases  forever  and  eventually  becomes  arbitrarily 
large.  Clearly  no  population  in  the  real  world  is  represented  by  this  model, 
since  no  limited  environment  can  sustain  an  infinite  population.   Thus  it  is 
necessary  to  adjust  the  model. 

Suppose  then  that  as  a  population  increases  above  the  maximum 
sustainable  capacity  of  the  environment,  the  competition  for  resources  causes 
some  of  the  population  to  die  off  or  be  killed.  Hence  a  particularly  large 
population  might  cause  the  population  to  actually  decrease  in  the  next  cycle. 
If  we  scale  our  model  so  that  0  <  x  <  1,  with  1  being  the  capacity  of  the 
environment,  then  incorporating  the  multiplicative  factor   (1  -  x)   into  our 
model  could  account  for  the  limitations  of  the  environment  on  large 
populations.  We  now  adjust  B  to  account  for  the  birthrate  (and  deathrate)  in 
a  fixed  environment  and  denote  its  new  value  by  X,  obtaining  the  refined 
model  xn+i  =  Xxn(l  -  xn).  This  model  shows  that  for  small  populations,  the 
population  growth  is  close  to  the  original  model  xn+i  =  Bxn.  However,  as  the 
population  approaches  the  maximum  sustainable  capacity  of  the 
environment,  its  growth  starts  to  decline  because  of  the  (1  -  xn)  factor. 

The  equation  xn+i  =  Xxn(l  -  xn),  X  >  0,  is  called  the  logistic  equation.  It  is 
perhaps  one  of  the  simplest  population  models  to  represent  the  interaction 
between  a  species  and  its  fixed  environment  (although  it  excludes 
interactions  with  other  species).  Allowing  x  to  vary  between  zero  and  one 
scales  the  model  to  any  population  size.  Restricting  the  parameter  0  <,  X  <  4, 
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we  see  below  that  the  model  never  allows  the  population  to  increase  above 
the  maximum  capacity  of  the  environment.   The  logistic  equation  has  been 
used  by  biologists  for  many  years,  and  it  has  been  a  major  goal  to  understand 
populations  that  are  modelled  by  this  equation.  One  methodology  has  been 
to  determine  experimentally  a  value  of  X   based  on  the  species  and  the 
environment,  and  then  to  observe  the  population  over  time  to  see  if  it 
behaves  according  to  the  model.  This  simple  model  is  a  fairly  good  predictor 
for  a  simple  species  like  bacteria  or  yeast  growing  in  a  culture,  but  it  has  not 
always  correctly  predicted  the  observed  results  for  more  highly  developed 
species,  like  mammals.  Some  populations  have  been  observed  to  settle  down 
to  a  fixed  value  for  all  time  increments,  others  have  been  seen  to  cycle 
between  two,  four,  or  even  larger  periods,  and  still  others  have  demonstrated 
completely  unpredictable  behavior  that  fit  no  known  model  at  all.  What 
happens  is  that,  for  certain  values  of  \,  the  logistic  equation  f(xn)  =  Xxn(l  -  xn) 
becomes  chaotic.  Thus  it  is  impossible  to  predict  its  behavior  even  after  a  few 
iterations.   We  next  examine  this  map  mathematically. 

First  note  that  for  values  of  xo  equal  to  0  or  1,  the  next  iteration  xi  =  0 
is  a  fixed  point  of  this  map,  since  f (0)  =  0.  Also,  for  values  of  xo  outside  the 
interval  [0, 1],  xi  becomes  negative  (which  does  not  make  physical  sense 
since  we  are  studying  populations).  Furthermore,  all  negative  values  diverge 
to  -o©  under  repeated  iteration  of  f.  Hence,  for  practical  reasons  we 
concentrate  on  the  dynamics  of  this  system  on  the  interval  [0, 1]. 

Suppose  0  <k  <  1.  Then  x  =  0  is  the  only  fixed  point  of  the  system  (see 
Figure  4.13).  Moreover,  all  values  of  xo  €  [0, 1]  tend  to  0  under  iteration  of 
f.  (Recall  that  if   I  f  (x)  I  <  1  for  a  fixed  point  x,  then  the  fixed  point  is 
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attracting.  Since  f  (x)  =  X(l  -  2x),  for  X  <  1,  we  obtain  f  (0)  =  X  <  1;  so  0  is  an 
attracting  fixed  point).  The  physical  interpretation  of  this  result  is  that  the 
population  does  not  reproduce  fast  enough  to  sustain  itself  and  eventually 
dies  out. 


f(x)  =  X 


►    X 


0  1 

Figure  4.13  The  logistic  equation  for  X  <  1. 


Now  increase  X  so  that  it  is  greater  than  1.  Then  0  becomes  a  repelling 
fixed  point  since  f  (0)  =  X  >  1.  However,  note  also  the  introduction  of  another 
fixed  point  p  in  the  interval  [0, 1],  given  by  p  =  1  -  1/X  (see  Figure  4.14).  If 
we  consider  values  of  X  such  that  1  <  X  <  3,  we  see  that  the  point  p  is  an 
attracting  fixed  point  since  I  f  (p)  I  =  I  X(l  -  2  +  2/X)  I  =  1 2  -  X I  <  1.  Hence,  for 
values  of  X  between  1  and  3,  0  is  a  repelling  fixed  point  and  p  =  l-l/A.  is 
an  attracting  fixed  point.  This  means  that  our  population  will  always  settle 
down  to  a  fixed  value  and  not  vary  in  cycles. 
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f(x)  =  x 


►    X 


Figure  4.14  The  logistic  equation  for  1  <  \  <  3. 


Now  let  A.  jump  to  the  value  X  =  4.  Since  xn+i  =  4xn(l  -  xn)  we  know 
that  this  map  is  chaotic  because  of  our  discussion  of  Newton's  method  for 
z2  =  -1.  We  now  demonstrate  this  idea  graphically  to  clarify  further  what  is 
actually  happening.  Consider  the  graph  of  f(x)  =  4x(l  -  x)  on  the  interval  [0, 
1]  in  Figure  4.15.  The  graph  represents  a  two-to-one  and  onto  function. 
Hence,  the  second  iterate  f2(x)  creates  a  condensed  copy  of  f(x)  on  each 
interval  [0,1/2]  and  [1/2,1].  Moreover,  that  f2(x)  has  four  fixed  points  (see 
Figure  4.16).  The  third  iteration  f3(x)  is  also  shown  in  Figure  4.16.  Observe 
that  we  are  creating  a  graphical  situation  similar  to  the  Baker  map  studied 
earlier  (with  parabolic  arcs  rather  than  straight  lines).  The  graphs  reveal  that 
f(x)   has  sensitive  dependence  on  initial  conditions,  topological  transitivity, 
and  a  dense  orbit;  so  it  is,  in  fact,  chaotic.  This  analysis  explains  why,  for  some 
populations,  the  observed  values  fluctuate  wildly,  despite  a  seemingly  simple 
mathematical  model  describing  the  behavior. 
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/  f(x)  =  x 


0  3/4     1 

Figure  4.15  The  logistic  equation  for  X  =  4. 


►  x 


►    X 


Figure  4.16  The  second  and  third  iterates  of  the  logistic  equation  for  X  =  4. 
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Another  way  to  show  that  the  map  f(x)  =  4x(l  -  x)  is  chaotic  is  to  use  the 
concept  of  topological  conjugacy.   Recalling  the  trigonometric  identity 
sin220  =  4sin29(l  -  sin28),  the  mapping  h(9)  =  sin20    shows  that    f(x)  =  4x(l  -  x) 
is  topologically  semi-conjugate  with  the  map  of  the  circle  g(9)  =  29.  This 
commutative  diagram  is  shown  in  Figure  4.17.    This  is  yet  another  way  to 
show  that  the  map  f(x)  =  4x(l  -  x)  is  chaotic. 


0 

—  g     -> 

29 

ih 

in 

sin  9 

—  f     -> 

sin220 

Figure  4.17  Topological  conjugacy  between  g(9)  =  29  and  f(x)  =  4x(l  -  x). 

Now  examine  what  happens  when  3  <  X  <  4.  We  observed  an  attracting 
fixed  point  for  X  <  3,  and  chaotic  behavior  for  X  -  4.  Thus  as  X  increases 
from  3  to  4  the  dynamics  of  the  system  must  change  dramatically  to 
produce  this  chaotic  behavior.  For  X  =  3,  the  magnitude  of  the  derivative  at 
the  fixed  point  p  =  1  -  1/X  =  2/3  is  equal  to  I  f  (2/3)  I  =  12  -  X\  =  1;  for  X  >  3, 
I  f  (p)  I  >  1.  This  means  that  as  the  value  of  X  passes  through  3,  the  fixed 
point  p  goes  from  attracting  to  repelling.  This  result  implies  that  both  of  the 
fixed  points,  0  and  p  =  1  -  1/X  are  repelling.  With  both  fixed  points 
repelling,  one  might  think  the  orbit  will  fluctuate  wildly  within  the  interval 
[0,1].  But  what  really  occurs  is  the  creation  of  an  attractive  two-cycle,  as 
shown  in  Figure  4.18.  This  splitting  of  a  fixed  point  into  an  attractive  cycle  is 
called  a  period-doubling  bifurcation. 
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Figure  4.18  An  attractive  two-cycle  in  the  logistic  equation  for  X  >  3. 

The  two-cycle  previously  described  is  one  of  the  phenomena  that  puzzled 
biologists  for  many  years.  Physically,  a  small  population  would  flourish  and 
breed  very  rapidly,  since  the  environment  could  easily  support  it.  After  a 
time  a  large  population  would  be  created  in  the  next  cycle.  This  population 
would  breed  more  sluggishly  because  it  was  nearer  the  capacity  of  the 
environment.    Furthermore,  the  population  would  not  approach  its  fixed 
value  (since  that  value  was  repelling).  Even  when  biologists  altered  the 
population  it  would  still  tend  back  to  the  attractive  cycle. 

As  we  continue  to  increase  the  value  of  X  beyond  X  =  3,  we  find  that  at 
X  ~  3.4495  the  two-cycle  becomes  repelling  and  an  attractive  four-cycle  is 
created.  Now,  as  with  X  -  3,  the  derivative  I  (f^VCp)  I   for  the  cyclic  points  pi 
and  p2  is  passing  through  1,  and  the  two-cycle  becomes  repelling;  another 
bifurcation  of  the  system  has  occured. 

If  we  continue  to  increase  X,  the  next  bifurcation  occurs  at  X  ~  3.544, 
where  the  four-cycle  becomes  repelling  and  an  attractive  eight-cycle  is  created. 
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Plotting  the  values  of  the  parameter  X  against  the  values  of  the  attractor  (the 
attracting  periodic  points  of  period  two,  four,  eight, . . .)  shows  how  each  cycle 
bifurcates  into  another  cycle  with  twice  the  period  as  the  previous  one.   Such 
a  bifurcation  diagram  is  shown  in  Figure  4.19  for  2.9  <  X  <  3.9.  The  parameter 
values  at  which  the  next  four  bifurcations  occur  are  X  =  3.5644,  3.5688,  3.5697, 
and  3.5699.  As  one  can  see,  these  numbers  are  getting  closer  together  and  they 
converge  to  a  value  where  the  map  becomes  chaotic. 


Figure  4.19  The  bifurcation  diagram  for  the  logistic  equation  as  X 
increases.  The  figure  is  from  Holden  (1986,  p.  46). 


The  apparent  convergence  of  these  numbers  led  to  the  discovery  of  a 
constant  that  appears  to  be  almost  universal  in  dynamical  systems.  Denote 
the  values  of  X  at  which  bifurcations  occur  by  Xk  for  k  =  1, 2, 3, ... ,  denote 
the  value  to  which  they  converge  by  Xoo,  and  consider  the  following.     Figure 
4.20  shows  the  second  iterate  of  the  logistic  equation.  Note  that  it  contains  a 
scaled  copy  of  the  first  iterate,  and  that  each  further  iterate  will  again  contain  a 
scaled  copy  of  the  original  map.  This  geometric  observation  provides  insight 
to  the  phenomenon  of  period-doubling.   At  X  =  Xk-i,  a  period-doubling 
occurs  in  the  system  that  lasts  for  a  duration  (in  X)  of  (Xk  -  Xk-l)-  The  next 
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period  doubling  occurs  at  X  =  Xk,  and  lasts  for  a  ^-duration  of  (Xk+1  -  Xk)-  It  is 
not  unreasonable  to  expect  a  "scaling"  in  the  successive  ^-lapses  based  on  the 
graphical  observation.  Hence,  if  we  assume  the  convergence  of  (Xoo  -  Xk)  to 
be  geometric,  then  A**  -  Xk  =  c/8k  where  c  is  a  constant  and  5  is  a  constant 
greater  than  one.  Using  this  to  solve  for  8  in  terms  of  the  Xk,  we  see  that 
8  =  (Xk  -  Xk-i)/(Xk+i  "  ^*  ^rom  our  previously  computed  values  of  Xk/  we 
see  that  8  =  4.6692,  which  is  called  the  Feigenbaum  constant,  named  after  the 
American  scientist  and  mathematician  Mitchell  Feigenbaum,  who  discovered 
it.    Using  this  constant  to  solve  for    Xoo,  we  find  that  Xco  *>  3.5699456.  Hence, 
with  A.  =  Xoo/  the  logistic  map  is  chaotic.  (As  an  aside,  the  Feigenbaum 
constant  has  been  found  in  many  different  dynamical  systems.  For  example, 
models  for  electrical  circuitry,  optical  systems,  and  economic  cycles,  as  well  as 
physical  systems  such  as  chemical  reactions,  erratic  heart  behavior,  and  even 
dripping  faucets,  have  all  exhibited  the  Feigenbaum  constant  when  viewed  in 
the  proper  phase  space.  With  its  value  now  known,  the  number  of 
dynamical  systems  where  the  Feigenbaum  constant  is  found  is  increasing 
quite  rapidly.) 


■►   x 


*►   x 


Figure  4.20  The  first  and  second  iterates  of  the  logistic  equation. 
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This  "period-doubling  route  to  chaos"  was  one  of  the  first  routes  to  chaos 
understood  by  mathematicians.   It  has  helped  biologists  understand  the  many 
different  results  associated  with  systems  modeled  with  the  logistic  map. 
While  this  discussion  about  bifurcations  and  the  period-doubling  route  to 
chaos  has  been  mostly  geometric  and  intuitive,  there  is  a  rigorous  underlying 
mathematical  theory  validating  these  concepts.   However,  this  theory 
requires  an  understanding  of  "kneading  theory/'  an  elaborate  version  of 
symbolic  dynamics  that  is  beyond  the  scope  of  this  thesis. 

While  the  logistic  map  is  very  restrictive  (in  that  it  models  only  a  single 
species  interacting  with  a  fixed  environment)  it  still  provides  a  great  deal  of 
information  about  population  dynamics.   When  the  model  is  further  refined 
with  the  introduction  of  another  species  (i.e.,  a  predator-prey  or  plant- 
herbivore  model),  the  possibility  for  chaotic  behavior  increases  with  the 
complexity  of  the  model.  Also,  the  logistic  equation  has  been  used  by  medical 
scientists  to  model  the  spread  of  infectious  diseases  (the  x  term  being  the 
number  who  are  contagious  and  the  (1  -  x)  term  being  the  number  who  have 
developed  immunity),  and  even  by  sociologists  to  model  the  spread  of 
rumors. 

We  now  continue  our  mathematical  analysis  of  the  logistic  map  beyond 
the  ranges  of  A.  that  correspond  to  physical  population  systems.  Specifically, 
we  examine  what  occurs  with  the  map  xn+i  =  focn(l  -  xn)  when  X  increases 
above  4.  The  simplest  way  to  observe  the  dynamics  is  graphically,  as  shown 
in  Figure  4.21.  As  before,  all  initial  values  outside  of  [0, 1]  diverge  to  -<*> 
under  infinite  iteration.  Notice  now,  however,  that  there  is  a  set  of  points 
inside  [0, 1]  mapped  outside  the  interval  after  the  first  iteration;  thus  they 
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diverge  to  •».  Examining  the  second  iteration,  we  see  that  the  two 
subintervals  remaining  after  the  first  iteration  are  mapped  to  [0, 1];  each  will 
then  contain  a  smaller  copy  of  the  original  map,  and  hence  a  subinterval 
mapped  outside  of    [0, 1]  (see  Figure  4.22).  In  fact,  if  we  consider  the  set  of 
points  which  do  not  eventually  diverge  to  -«>,  by  this  construction  we  obtain 
a  Cantor  set  of  points  remaining  in  the  interval  after  infinite  iteration  for 
X  >  4.  Values  of  the  interval  [0, 1]  mapped  outside  of  the  unit  interval,  and 
hence  diverging  to  -«>,  are  said  to  be  values  that  escape  under  iteration  of  f. 
The  value  of  X  for  which  the  interval  (1/3,2/3)  escapes  under  the  first 
iteration  of  f  can  be  found  by  letting  1  ■  X(l/3)(1  -  1/3),  which  yields  X  =  9/2. 
However,  iterating  the  logistic  equation  with  X  =  9/2  will  not  yield  the 
classical  Cantor  set,  since  the  subintervals  created  in  the  second  iteration  will 
be  skewed  and  will  not  be  exact  middle-thirds  intervals.  In  fact,  the  point 


f(x) 


0  escape 

Figure  4.21  The  logistic  equation  with  X  >  4  and  the  interval  that  escapes. 


107 


KT^ 


Figure  4.22  The  second  iterate  of  the  logistic  equation  with  X  >  4. 

x  =  1/9  escapes  under  the  second  iteration.  However,  we  can  also  see  how  a 
Cantor-like  set  of  points  remains  in  the  set.  Considering  X  -  9/2,  we  see  that 
f(l/3)  =  1,  and  since  f(l)  -  0,  with  0  being  a  fixed  point  of  the  map,  the  point 
1/3  forever  remains  fixed  in  the  set.  For  X  >  4,  the  logistic  equation  is  chaotic 
on  a  Cantor  subset  of  [0, 1],  while  on  the  complement  of  this  set,  all  orbits 
diverge  predictably.  This  contrasts  with  the  logistic  equation  with  1  =  4, 
which  is  chaotic  on  the  entire  interval  [0,1]. 

L     BIFURCATIONS 

In  the  last  section  we  introduced  the  concept  of  bifurcations  as  period- 
doubling  phenomenon  leading  to  chaos  in  the  logistic  map.  We  now  study 
the  general  concept  of  bifurcation  and  its  interpretation  in  the  real  world. 
The  mathematical  theory  in  this  chapter  is  from  Rasband  (1990), 
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Guckenheimer  (1990),  and  Seydel  (1988),  while  the  particular  definitions  and 
theorems  are  from  Devaney  (1989).  The  examples  are  taken  from  Moon 
(1987)  and  Holden  (1986). 

Bifurcations  are  normally  associated  with  a  change  in  the  physical  state  of 
a  system.  In  population  dynamics,  one  bifurcation  is  exhibited  by  a 
population  changing  from  tending  towards  an  attracting  fixed  point  to 
tending  towards  an  attractive  cycle.  Changes  in  the  states  of  other  physical 
systems  (such  as  water  cooling  until  it  freezes,  a  beam  bending  under  a  load 
until  it  buckles,  a  balloon  being  blown  up  until  it  bursts,  or  even  the  crash  of 
the  stock  market)  are  also  examples  of  bifurcations.  Clearly,  some  of  these 
behaviors  have  such  complex  and  intricate  mathematical  models  that  it 
would  be  virtually  impossible  to  quantify  them  and  determine  the  exact 
nature  of  the  bifurcation.  Nevertheless,  there  are  several  systems  which  have 
been  modeled  successfully  and  for  which  the  mathematical  analysis  of 
bifurcation  accurately  predicts  gross  changes  in  the  state  of  the  system. 

Mathematically,  there  are  several  types  of  bifurcations.  First  we  discuss 
those  that  occur  in  one  dimension  (as  with  our  previous  example  of  the 
logistic  equation).  In  the  logistic  equation,  the  bifurcation  occurs  as  the 
parameter  value  X  is  varied  until  it  passes  through  certain  values  which 
changed  dramatically  the  behavior  of  the  system.  This  example  suggests  we 
consider  families  of  functions  of  real  variables,  where  the  functions  depend 
"smoothly"  on  a  parameter.  Thus,  define  g(x,  X)  =  fx(x),  where  f(x)  is  a  C°° 
function  of  x  for  each  fixed  X,  and  g(x,  X)  depends  smoothly  on  X.  The 
logistic  equation  g(x,  X)  =  fx(x)  =  Xx(l  -  x)    is  one  such  family. 
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We  first  present  a  general  result  about  bifurcation  theory  that  applies  to 
all  one-dimensional  bifurcations.   Then  we  give  specific  examples  of  typical 
bifurcations  and  physical  systems  that  behave  similarly  in  the  real  world.  The 
general  result  is  that  bifurcations  occur  only  near  non-hyperbolic  fixed  and 
periodic  points.  (Recall  that  a  hyperbolic  periodic  point  is  one  where 
I  (fn)'(p)  I  *  I,  where  n  is  the  prime  period  of  the  point;  for  a  fixed  point, 
n  =  1.)  While  the  theory  applies  to  both  fixed  and  periodic  points,  we  present 
it  here  only  for  fixed  points.  Recalling  further  that  a  fixed  point  p  is  an 
attractor  if   I  f  (p)  I  <  1,  and  a  repellor  if  I  f  (p)  I  >  1,  we  would  expect  to  find 
interesting  behavior  as    I  f  (p)  I   passes  through  1.  Such  points  p  are  the  only 
points  where  bifurcations  can  and  do  occur. 

BIFURCATION  THEOREM.  Let  fx  be  a  one  parameter  family  of  functions, 
with  p  a  hyperbolic  fixed  point  for  some  fixed   Xq,  i.e.,   fxo(p)  ■  P   and 
f  Xo(p)  *  1-  Then  there  exist  intervals,  I  about  p  and  N  about  Ao,  and  a 
smooth  function  h:  N->I  such  that  h&o)  =  p  and  fx(h(A.))  =  h(X).  Moreover, 
f\  has  no  other  fixed  points  in  I. 

PROOF.  Consider  the  function  defined  by   G(x,  X)  =  fx(x)  -  x.  By  hypothesis, 
G(p,  Xo)  =  0,  and  dG/dx  I  (p,  Xo)  =  *Xo'(p)  -1*0.  By  the  Implicit  Function 
Theorem,  there  are  intervals,  I  about  p  and  N  about  Xo,  and  a  smooth 
function  h:  N->I  such  that  hfto)  =  p,  and  G(h(X.),  X)  =  0  for  all  X  e  N. 
Moreover,  g(x,  X)  *  0  unless  x  =  h(X).  This  completes  the  proof. 
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A  graphical  demonstration  of  this  theorem  is  shown  in  Figure  4.23.   The 
graph  shows  that  nearby  graphs  fxi  and  fx2  must  have  the  same  property  for 
a  sufficiently  small  interval  about  Xo,  because  the  graph  of  i\  meets  the  line 
y  =  x  at  an  angle  at  (p,  p)  and  since  f  varies  smoothly  with  X.  Hence,  there  is 
only  one  fixed  point  near  p  for  X  in  some  neighborhood  of  Xo-  At  this 
point  there  is  a  function  g\  that  is  topologically  conjugate  to  i\  via  the  map 
Hx(x)  =  x  -  h(A.)  for  which  the  origin  is  always  a  fixed  point.  This  fact  allows 
us  to  study  maps  with  0  as  a  fixed  point,  and  apply  the  results  to  any  map  i\ 
with  a  nonzero  fixed  point.  (To  see  this,  consider  fx  with    f(h(X))  =  i(X).     If 
gx(x)  =  fx(x  +  UX))  -  h(X),  then  gx(0)  =  (x(h(X))  -  HX)  =  0  for  all  X,  so  0  is 
always  a  fixed  point  for  g.)  Thus  it  suffices  to  present  our  complete  results  for 
one-dimensional  bifurcations  using  functions  whose  fixed  point  is  the 
origin. 


►  x 


Figure  4.23  Schematic  of  the  Bifurcation  theorem. 
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Now  that  we  know  that  bifurcations  occur  at  non-hyperbolic  periodic 
points,  we  are  prepared  to  study  specific  types  of  one  dimensional 
bifurcations.  The  first,  and  perhaps  mathematically  the  simplest  type,  is  the 
tangent  bifurcation  (also  called  the  saddle-node).  The  family  fx(x)  =  Xex,  for 
X  >  0,  has  a  bifurcation  at  x  =  0  when  X  =  1/e.  If  X  >  1/e  the  function  f\  has 
no  fixed  points,  and  f\n(x)  ->  ©°  for  all  x.  When  X  =  1/e,  the  function  has  one 
fixed  point  at  p  =  1.  This  point  is  attracting  from  the  left  and  repelling  from 
the  right.  This  last  condition,  which  may  occur  for  a  non-hyperbolic  fixed 
point,  is  called  semi-stable.  For  values  of  X  <  1/e,  i\  has  two  fixed  points,  pi 
and  p2,  with  pi  attracting  and  p2  repelling.  Graphs  of  these  three  cases  are 
displayed  in  Figure  4.24,  and  the  graphical  analyses  are  left  to  the  reader. 


Figure  4.24  Graphs  of  f(x)  =  Xe*  with  (1)  X  <  1/e  (2)  X  =  1/e  (3)  X  >  1/e. 

We  can  see  that  this  bifurcation  occurs  as  the  graph  of  fx  becomes 
tangent  to  the  line  y  =  x.  First  one,  then  two,  fixed  points  emerge  as  the 
parameter  passes  through  the  critical  value.  We  can  also  verify  that  this 
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bifurcation  occurrs  at  a  non-hyperbolic  fixed  point  p  =  1,  since    fi/e'(l)  ■ 
(l/eje1  =  1.  This  bifurcation  gets  its  name  from  the  way  fx  approaches  the 
line  y  =  x  tangentially.  The  bifurcation  diagram  for  the  function  i\/e, 
plotting  the  fixed  points  p  on  the  vertical  axis  versus  X  is  shown  in  Figure 
4.25. 


1/e  '(not  to  scale) 

Figure  4.25  Bifurcation  diagram  for  f (x)  =  Xe*. 


The  mathematical  conditions  guaranteeing  the  occurrence  of  a  tangent 
bifurcation,  and  the  resulting  bifurcation,  are  as  follows. 

THE  TANGENT  BIFURCATION.  Given  that 

1.  ftu>(0)  =  0, 

2.  fxo'(O)  =  1, 

3.  fxo"(0)  *  0,  and 

4.  afx/9Mx=Xo*0. 
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Then  there  exists  an  interval  I  about  0,  and  a  smooth  function  h:  I— »R 
satisfying  h(0)  =  Xof  such  that  fh(x)(x)  =  x.  Moreover,  h'(0)  =  0  and  h"(0)  *  0. 

The  tangent  bifurcation  has  received  considerable  attention  of  engineers 
recently  because  it  is  known  that  a  beam  or  small  arch  bending  under  a  load 
undergoes  a  tangent  bifurcation  at  its  buckling  point.  Additionally, 
astronomers  have  identified  the  equilibrium  response  of  massive  cold  stars 
with  the  tangent  bifurcation. 

A  bifurcation  familiar  to  us  from  the  analysis  of  the  logistic  equation  is 
the  period-doubling  bifurcation.  We  now  give  another  example  of  the 
period-doubling  bifurcation.  Consider  the  family  of  functions  f\(x)  =  Xex,  for 
X  <  0.  Graphs  of  these  functions  are  shown  in  Figure  4.26.  When  X  =  -e, 
fx(-l)  =  -1,  and   I  fx'(-l)  I  =   I  -eCe"1)  1  =  1-11=1;  thus  p  =  -1  is  a  non-hyperbolic 
fixed  point.  We  would  expect  a  bifurcation  to  occur  as  X  passes  through  the 
critical  value  -e,  and  indeed  one  does  occur.  When  X  >  -e,  the  fixed  point  p 
is  attracting,  and  when  X  <  -e,  it  is  repelling.  Hence,  the  nature  of  the  fixed 
point  changes  as  X  passes  through  -e,  but  this  is  not  all  that  occurs.  When  X 
<  -e,  the  graph  of  fx2(x)  is  an  increasing  function  that  is  concave  up  if 
fx(x)  <  -1,  and  concave  down  if  f^(x)  >  1.  Since  fx2(x)  has  two  fixed  points, 
this  corresponds  to  a  two-cycle  in  fx(x).  Figure  4.27  shows  the  bifurcation 
diagram  for  the  above  example,  again  plotting  the  fixed  point  and  the  periodic 
points  on  the  y-axis  against  the  parameter  X. 
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Figure  4.26  Graphs  of  the  family  fx(x)  =  Xex  for  X  >  -e  and  X  <  -e. 


fixed  point 
is    attracting 


fixed   point 
is  repelling 


two-cycle 
Figure  4.27  Bifurcation  diagram  for  f(x)  =  Xe*  with  X  <  0. 

In  a  period-doubling  bifurcation,  as  the  parameter  passes  through  its 
critical  value,  the  attracting  fixed  point  becomes  repelling  and  a  cycle  of  period 
two  emerges.  Of  course,  as  we  saw  with  the  logistic  equation,  a  bifurcation 
can  also  occur  that  makes  an  attracting    n-cycle  repelling  and  creates  an 
attracting    2n-cycle;  hence  the  name  "period-doubling."  The  formal 
mathematical  characterization  of  the  period-doubling  bifurcation  is  given  in 
the  following  result. 
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THE  PERIOD-DOUBLING  BIFURCATION.  Given  that 

1.  fx(0)  =  0  for  all  X  in  an  interval  about  ^o, 

2.  fxo'(O)  =  -1,  and 

3.  d(tf)7dMuXo(0)*0. 

Then  there  is  an  interval  I  about  0,  and  a  function  h:  I-»R,  such  that 
fh(x)(x)  *  x,  but  fh(x)2(x)  -  x. 

As  shown  earlier,  without  loss  of  generality  we  are  able  to  consider  functions 
with  fixed  points  at  the  origin  and  apply  the  same  results  in  the  general  case 
for  arbitrary  fixed  points.  The  proofs  of  both  results  concerning  bifurcations 
can  be  found  in  Devaney  (1989),  and  involve  little  more  than  application  of 
the  Implicit  Function  Theorem  and  knowledge  of  partial  derivatives. 

A  situation  where  period-doubling  bifurcations  have  been  observed 
occurs  in  the  field  of  electrocardiology.  Electrochemical  events  in  the  heart, 
monitored  by  electrocardiograms  (ECG),  show  periodic  activity  within  the 
atria  and  ventricles  of  the  heart.  Abnormal  cardiac  rhythms,  such  as 
arrhythmia,  have  long  been  referred  to  as  "chaotic  heart  action"  (in  the 
descriptive,  non-mathematical  sense)  by  cardiologists.  Mostly,  these 
fibrillations,  or  chaotic  heart  rhythms,  have  been  recorded  on  ECGs 
immediately  prior  to  a  patient's  death.   By  concentrating  on  the  heart's 
natural  pacemakers,  called  ectopic  foci  and  which  are  located  throughout  the 
heart,  cardiologists  have  found  that  in  arrhythmia,  the  rhythm  of  the  heart 
undergoes  a  series  of  period-doubling  bifurcations  leading  eventually  to 
chaotic  behavior  (in  the  precise  mathematical  sense).   While  the  human 
heart  requires  a  much  more  complex  model  than  the  one-dimensional 
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dynamical  system,  there  is  a  high  correlation  between  the  period-doubling  we 
have  studied  and  the  period-doubling  on  the  ECG.  Cardiologists  are  using 
this  new  theory  to  try  to  explain  what  actually  happens  in  the  heart  in  order 
to  prevent  or  cure  the  condition.   A  current  solution  to  this  problem  is  to 
identify  patients  who  are  susceptible  to  arrhythmia  and  implant  in  their 
chests  a  device  that  detects  fibrillations  and  gives  the  heart  an  electronic 
"kick"  out  of  its  period-doubling  path.  This  device  has  proven  successful,  but 
it  has  been  difficult  to  identify  the  patients  who  will  benefit  from  its  use. 

Turbulence  is  one  of  the  classical  and  unexplained  phenomena  in 
physics.  In  the  past  it  has  been  so  inaccessible  to  physicists  and  engineers  that 
systems  have  been  designed  normally  with  "fudge"  factors,  or  factors  of  safety, 
designed  to  compensate  for  the  effects  of  turbulence.  However,  the  issue  of 
turbulence  has  rarely  been  directly  taken  on.  Recently,  in  a  physical  model 
that  creates  turbulence  in  a  very  simple  fluid  setting,  scientists  have  observed 
period-doubling  bifurcations  in  the  fluid  leading  directly  to  the  onset  of 
turbulence.  These  observations  may  be  the  first  steps  towards  an 
understanding  of  this  elusive  behavior. 

A  few  more  examples  of  how  a  changing  parameter  in  the  real  world  can 
lead  to  bifurcations  will  emphasize  the  importance  of  studying  this  concept. 
In  the  logistic  equation  it  is  reasonable  to  assume  that  the  parameter  X,  which 
measures  the  species'  interaction  with  its  environment,  does  not  change  too 
drastically  from  cycle  to  cycle.  On  the  other  hand,  in  a  fluid  system,  it  is  just 
as  reasonable  to  see  how  the  changing  velocity  of  a  flow  can  change  a 
parameter  quickly  enough  to  induce  the  period-doubling  bifurcations  that 
lead  to  turbulence.  So  it  might  seem  reasonable  to  surmise  that  the  panel 
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flutter  on  a  supersonic  aircraft  could  be  avoided  through  a  careful  design  that 
stays  away  from  dangerous  parameter  values.  However,  not  every 
atmospheric  condition  can  be  predicted  accurately  by  a  model,  or  duplicated  in 
a  wind  tunnel.   Moreover,  atmospheric  phenomena,  such  as  ice  buildup  on 
wings  or  wind  shears,  can  quickly  push  parameters  into  critical  regions.  In 
complex  systems  one  can  rarely  anticipate  every  set  of  parameter  values, 
which  is  an  important  reason  to  pay  careful  attention  to  parameter  space. 

J.     SARKOVSKIFS  THEOREM 

We  now  present  another  remarkable  result  about  one-dimensional 
dynamical  systems  due  to  the  Russian  mathematician,  A.  N.  Sarkovskii.   The 
development  here  is  from  the  article  by  Devaney  (1989).  First,  order  the 
positive  integers  in  the  following  manner: 

3»5»7»9»'"» 

2-3  »  2-5  »  2-7  »  2-9  » •  •  •  » 

22-3  »  22-5  »  22-7  »  22-9  »  •  • » 

23-3  »  23-5  »  23-7  »  23-9  » •  • » 

»  •  •  • »  2r  »  23  »  22  »  2  »  1. 
Here,  we  have  first  listed  the  odd  numbers  in  ascending  order,  followed  by 
two  times  the  odd  numbers,  then  22  times  the  odd  numbers  and  so  on, 
through  every  positive  integer  power  of  2  times  the  odd  numbers,  finally 
followed  by  the  powers  of  2  in  descending  order.  Given  this  ordering  of  the 
positive  integers,  we  present  Sarkovskii's  theorem. 
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SARKOVSKII'S  THEOREM.   Suppose  f:  R-»R  is  continuous.  If  f  has  a 
periodic  point  of  period  n,  then  f  also  has  a  periodic  point  of  period  k  for 
all  k  with  n  »  k  in  the  above  ordering. 

While  this  result  only  holds  for  one-dimensional  systems,  it  is 
remarkably  powerful  due  to  its  lack  of  hypothesis.  In  fact,  the  only 
requirement  on  f  is  that  it  be  continuous. 

Rather  than  proving  Sarkovskii's  theorem  itself,  we  prove  a  special  case 
in  what  follows.  This  proof  is  found  in  the  article  by  Devaney  (1989)  and 
differs  from  the  complete  proof  of  Sarkovskii's  theorem  by  requiring  less 
bookkeeping.  The  full  proof  of  Sarkovskii's  theorem  can  be  found  in 
Devaney  (1989). 

PROPOSITION.  Suppose  a  continuous  map  f:  R-»R  has  a  cycle  of  period  3. 
Then  f  has  periodic  points  of  all  periods. 

PROOF.  The  proof  depends  on  two  observations.  First,  if  I  and  J  are  closed 
intervals  with  J  z>  I  and  f (D  3  J,  then  f  has  a  fixed  point  on  I.  (This  is 
similar  to  the  result  about  fixed  points  we  proved  in  Section  EC).  The  second 
observation  is  as  follows:  suppose  Aq,  Ai,  . . . ,  An  are  closed  intervals  such 
that  f(Ai)  3  Aj+i  for  i  =  0, 1, . . . ,  n-1.  Then  there  exists  at  least  one 
subinterval  Jo  of  Aq  which  is  mapped  onto  Ai.  Similarly,  there  is  a 
subinterval  Ji  of  Ai  mapped  onto  A2,  and  hence  a  subinterval  Ji  of  Jo  such 
that  Ai  z>  f(Ji)  and  f2(Ji)  =  A2.  Continuing  in  this  manner,  we  obtain  a 
nested  sequence  of  subintervals  which  map  into  the  various  Ai ,  in  order. 
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To  prove  the  proposition,  let  a,  b,  and  c  be  such  that  f(a)  =  b,  f(b)  =  c,  and 
f(c)  =  a.  With  only  a  slight  loss  of  generality,  assume  a  >  b  >  c  (the  other  case 
is  handled  similarly).   This  hypothesis  is  shown  in  Figure  4.28. 


f(x) 


Z ►  x 

c  b        a 

Figure  4.28  A  three-cycle  with  f(a)  ■  b,  f (b)  ■  c,  f  (c)  =  a,  and  a  >  b  >  c. 


We  let  Io  -  [b,  a]  and  Ii  =  [c,  b].  By  our  assumptions,  f(Io)  z>  Ij  and 
ffli)  3  Io  u  Ii.  Figure  4.27  shows  there  is  a  fixed  point  between  c  and  b. 
Similarly,  f2  has  fixed  points  between  a  and  b,  so  at  least  one  of  these  points 
must  be  of  period  2.  Fixing   n  >  3,  we  now  produce  a  cycle  of  period  n. 

We  first  find  a  nested  sequence  Ao,  Ai, . . . ,  An_2,  of  subintervals  of  Ii  as 
follows:  let  Aq  =  Ii-  Since  ffli)  3  Ii,  there  is  a  subinterval  Ai  of  Ao  such 
that  f(Ai)  =  Ao  =  Ii.      By  induction,  we  can  find  a  subinterval  An_2  of  An_3 
such  that  f(An_2)  =  An_3,   f2(An-2)  =  An^, . . . ,  and  fr-HAad  ■  Ao  =  Ii.  Since 
f(Il)  3  Io,  there  is  a  subinterval  An-i  of  An-2  such  that  P^CAn-i)  =  Io-  Finally, 
since  f(Io)  3  Ii,  we  have  PKAn-i)  3  Ii  3  An_i.  Hence,  f*1  has  a  fixed  point  in 
An_i  from  the  first  observation.  We  now  show  this  point  has  period  n  under 
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f.  Since  Ii  =>  ^(An-i)  for  i  =  0, 1, . . . ,  n-2,  but  Io  z>  fn-KAn.i),  and    f*KAn)  =>  l\, 
this  point  has  its  first  n-2  iterates  in  Ii,  then  jumps  to  Io,  in  the  n  -  1 
iteration,  and  finally  back  to  1\.  This  completes  the  proof. 

In  addition  to  the  above  proposition,  another  corollary  to  Sarkovskii's 
theorem  states:  //  /  has  a  periodic  point  which  is  not  a  power  of  2,  then  f 
has  infinitely  many  periodic  points. 

Sarkovskii's  theorem  provides  considerable  information  about  a 
function.    For  example,  it  would  be  very  difficult  to  check  directly  whether 
the  function  f(x)  =  1  +  (5/2)x  -  (3/2)x2  is  chaotic.  However,  since  f(0)  =  1, 
f(l)  =  2,  and  f(2)  =  0,  the  function  has  a  three-cycle.  Thus  Sarkovskii's 
theorem  tells  us  it  has  cycles  of  all  periods.  Hence,  we  automatically  know  it 
is  chaotic  without  having  to  check  the  three  defining  conditions. 

If  we  consider  the  period-doubling  bifurcations  as  a  route  to  chaos,  then 
only  finitely  many  periodic  points  must  have  the  periods  1,  2,  22, 23, 24, . . . , 
2N  for  some  N.  Then  as  the  parameter  varies  and  the  dynamics  of  the 
system  become  more  complex,  we  introduce  periods  in  a  specific  order  2N+1, 
2N+2, ....  This  argument  does  not  claim  that  the  new  orbits  appear  as  period- 
doubling  bifurcations,  but  that  something  similar  must  occur. 

We  cannot  derive  a  converse  to  the  theorem  from  the  special  ordering.  If 
we  find  a  cycle  of  period  k,  and  n  »  k  in  the  ordering,  there  is  no  guarantee 
that  a  cycle  of  period  n  exists. 

While  we  used  a  graphical  analysis  to  show  that  the  Baker  map 
f(x)  =  2x  (modi)   is  chaotic,  we  now  confirm  this  using  Sarkovskii's  theorem. 
Although  the  Baker  map  is  not  continuous,  since  every  iterate  has  fixed 
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points  on  [0, 1],  we  can  apply  the  above  proof  of  the  proposition  directly  to  it. 
Thus,  if  we  can  find  a  three-cycle  of  the  Baker  map,  we  will  know  it  is  chaotic. 
But  that  is  easy:  since  the  points  1/7, 2/7,  and  4/7  form  such  a  cycle,  the 
Baker  map  is  chaotic. 

Another  example  using  Sarkovskii's  theorem  to  find  chaotic  behavior  is 
with  the  function  f(x)  =  x2  +  c,  where  c  =  -1.755.  Using  a  computer,  we  can 
verify  the  attracting  orbit  of  period  3  given  by  f(0)  =  -1.755,  f (-1.755)  * 

1.325  (to  four  decimal  places),  and  f (1.325)  =  0.  This  three-cycle  guarantees 
cycles  of  all  periods  and  an  infinite  number  of  periodic  points.  Nevertheless, 
regardless  of  the  initial  value,  the  orbit  is  always  attracted  to  this  three-cycle. 
But  then  where  are  the  other  periodic  points  if  only  three  are  found  by 
computational  iteration?  The  answer  is  that  all  other  cycles  are  repelling  and, 
because  of  computer  round-off  error,  iterates  always  move  away  from  a 
repelling  cycle  (unless  the  points  are  rational  numbers  represented  exactly  up 
to  the  precision  of  the  computer).  This  example  further  illustrates  the 
sensitive  dependence  on  initial  conditions  of  chaotic  dynamical  systems. 

Now  why,  of  all  the  infinite  orbits,  is  only  one  of  them  attracting?  This 
question  is  extremely  complicated  in  general.  However,  in  the  case  of  our 
specific  example  it  can  be  answered  using  a  result  from  complex  analysis.  For 
a  complex  analytic  map,  each  attracting  orbit  attracts  at  least  one  critical  point 
for  the  map.  Now  the  map  f(z)  =  z2  -  1.755,  is  an  analytic  function  with  only 
one  critical  point  z  =  0.  Hence  the  function  has  only  one  attracting  orbit. 
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K.    THE  QUADRATIC  FAMILY  REVISITED 

The  last  example  showed  another  member  of  the  quadratic  family  of 
maps;  namely,  f(x)  =  x2  +  c.  We  now  study  this  map  in  greater  detail  as  we 
vary  the  parameter  c.  The  results  in  this  section  are  from  Falconer  (1990)  and 
the  article  by  Devaney  (1989). 

First  note  that  for   c>  1/4  the  graph  of  f(x)  =  x2  +  c  lies  above  the 
diagonal  y  =  x  and  PXx)-**  for  all  xe  R.  For  c  =  1/4,  the  graph  of  f(x)  is 
tangent  to  the  diagonal  when  x  =  1/2    (which  is  a  single  fixed  point).  Finally, 
for  c<l/4,  f  has  two  fixed  points  which  we  denote  by  pi  and  p2   where 
pi  <  p2-  This  is  an  example  of  the  tangent  bifurcation  as  c  passes  through 
the  value  1/4.  These  three  cases  are  shown  in  Figure  4.29. 


►  x 


Figure  4.29  Graphs  of  f(x)  =  x2  +  c  for  c>  1  /4,  c  =  1  /4,  and  c  <  1  /4. 
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Now  observe  that  for  all  c<l/4,  p2>l/2/so  If(p2)l  =  1 2p2 1  >  1- 
However,  for  -3/4  <  c  <  1/4,  I  f  (pi)  I  <  1  since  -1/2  <  pi  <  1/2.  Finally  for 
c  <  -3/4,    I  f'(pi)  I  >  1  since  pi  <  -1/2.  This  demonstrates  that  after  the  tangent 
bifurcation  occurs  (as  c  passes  through  1/4)  the  fixed  point  pi  is  attracting 
and  p2  is  repelling.  However,  as  c  passes  through  -3/4,  pi  also  becomes 
repelling  (and  in  fact,  a  two-cycle  forms  around  pi).  The  graphs  of  f(x)  =  x2  +  c 
as  c  passes  through  -3/4  are  shown  in  Figure  4.30.  As  c  continues  to 
decrease,  we  get  another  sequence  of  period-doubling  bifurcations,  similar  to 
those  we  saw  with  the  logistic  equation.  The  frequency  of  these  bifurcations  is 
also  governed  by  the  Feigenbaum  constant. 


►  x 


C  >  -3/4  / 


C  <  -3/4 


Figure  4.30  The  graph  of  f(x)  =  x2  +  c  as  c  passes  through  -3/4. 

We  have  developed  three  ways  to  determine  if  this  period-doubling  leads 
to  chaos.  The  first  is  to  graph  f(x),  and  f^x)  through  f"(x),  for  certain  values 
of  c  and  observe  the  fixed  points  of  these  graphs.  The  second  way  is  to  solve 
algebraically  for  the  roots  of  f^x)  =  x  and  observe  the  fixed  points  as  the  roots 


124 


of  these  polynomials.  For  the  case  f2,  solving  the  equation  f2(x)  =  x  yields  the 
fourth-order  equation  (x2  +  c)2  +  c  =  x,  or  x4  +  2cx2  -  x  +  c2  +  c  =  0.  Since  we 
already  know  two  of  the  roots  (pi  and  P2)  for  any  given  value  of  c,  we  can 
solve  directly  for  the  other  two.  However,  as  n  increases,  it  becomes  very 
difficult  to  solve  for  the  roots  of  PKx).  Third,  we  can  simply  find  a  three-cycle 
for  certain  values  of  c  (for  example  c  =  -1.755)  and  appeal  to  Sarkovskii's 
theorem  as  proof  that  periods  of  all  other  orders  do  exist. 

Continuing  to  analyze  the  family  f(x)  =  x2  +  c,  we  see  from  Figure  4.28 
that  for  all  c  <  1  /4,  if  I  xo  I  >  p2,  then  f^Kxo)-*00.  Thus,  we  can  focus  our 
attention  on  the  interval  I  =  [-p2,  P2]  where  all  of  the  interesting  dynamics 
occur.  Let  us  further  restrict  our  attention  to  the  range  of  parameter  values  c 
<  -2.  If  we  analyze  this  function  on  a  computer  for  almost  any  initial  value, 
the  iterates  of  f  to  go  to  infinity.  However,  as  shown  below,  there  are  many 
orbits  which  do  not  escape  under  iteration  of  f. 

The  graph  of  f  for  c  <  -2  is  shown  in  Figure  4.31.  If  we  consider  the 
interval  I  =  [-p2,  p2l  note  that  there  is  a  subinterval  Ao  of  I  that  maps  to 
values  outside  of  I,  hence  all  points  in  Ao  will  escape  to  <*>.  Consider  the 
graph  of  f2,  a  similar  analysis  to  the  one  performed  for  the  logistic  map  with 
X  >  4  shows  that  two  subintervals  of  I  -  Ao  are  again  mapped  outside  of  I,  so 
again  points  escape  to  oo.  The  second  iteration  of  f  on  the  interval  I  is 
shown  in  Figure  4.32.  The  subintervals  that  escape  I  on  the  second  iteration 
are  labeled  Ai  and  A2  in  the  figure. 

Analyzing  the  points  that  remain  in  I    after  infinite  iteration,  we  deduce 
a  Cantor  set  has  been  constructed  on  I.  While  it  is  not  necessarily  the  classical 
Cantor  set,  nevertheless,  it  contains  no  intervals  despite  an  infinite  number 


125 


(P2.  P2) 


►    X 


Figure  4.31  The  graph  of  f  (x)  =  x2  +  c  for  c  <  -2. 


(P2.   P2) 


Figure  4.32  The  graph  of  f2(x)  for  c  <  -2. 
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of  points.  (As  we  mentioned  in  Section  IIE,  this  could  be  an  example  of  a  "fat" 
Cantor  set,  depending  on  the  value  of  c).  We  now  define  the  set 

A  =  {x e  I:  fn(x)  dVn^O), 
and  assert  that  A  is  a  Cantor  set. 

The  dynamics  of  f  on  R  -  A  are  quite  simple  because  every  initial  value 
tends  to  +<*>  under  infinite  iteration.  We  want  to  know  what  happens  to  the 
set  A.  To  determine  this,  we  simplify  the  analysis  through  symbolic 
dynamics.  Recall  that  the  shift  map  on  code  space  a:  X-»X  is  a  continuous 
mapping.  We  now  try  to  relate  a  and  f.  When  we  remove  the  interval  Ao 
from  I,  two  subintervals  remain,  denoted  by  Io  and  Ii  (see  Figure  4.33). 
Hence,  if   x  e  A,  then  f"(x)  e  Io  u  Ii    for  all  n  £  0.  Next  define  the  itinerary  of 
x  by  S(x)  =  (sosiS2.  •  •)  where  Si  e  {0, 1}  and  Si  =  k  if  and  only  if  fj(x)  e  1^. 


Figure  4.33  The  subintervals  Io  and  l\. 
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We  now  show  that  S:  A-»X  is  a  homeomorphism.  First,  to  see  that  S  is 
one-to-one,  let  x,  y  e  A  and  suppose  S(x)  =  Sty).  Then  for  each  n,  fn(x)  and 
fn(y)  both  lie  on  the  same  side  of  Ao-  It  follows  that  f  is  monotonic  on  the 
interval  between  fn(x)  and  f^Ky).  Hence,  all  points  in  this  interval  remain  in 
Io  u  Ii.  This  observation  contradicts  the  fact  that  A  is  a  Cantor  set  and 
contains  no  intervals,  so  S  is  one-to  one. 

To  see  that  S  is  onto,  for  a  closed  interval  J,  set  f_n(J)  =  {x  e  I:  fn(x)  e  J).  If 
I  z>  J,  then  f*1(J)  consists  of  two  subintervals,  one  in  Io  and  one  in  Ij.  Let 
S€  I  with  s  =  S0S1S2.  • . ,  and  define 

Iso...sn  -  (x  €  I:  x  €  Iso,  f(x)  e  L;i, . . . ,  f*(x)  e  Isn} , 
so  Iso...sn  ■  Iso  ri  f_1(Isi)  n . . .  n  f-n(Isn).  We  claim  that  the  Iso...sn  form  a 
nested  sequence  of  nonempty  closed  intervals  as  n-»<*>.  Note  that 
Iso...sn  =  Iso  n  f"1Ctsl...sn)«  By  induction,  we  assume  that  Isi...sn    is  a  nonempty 
closed  subinterval  so  that  f^Gsi^sn)  consists  of  two  subintervals,  one  in  Io 
and  one  in  l\.  Hence,  ISo...sn  is  a  single  closed  subinterval.  These  intervals 
are  nested  since  ISo...sn  =  Iso...sn-l  n  f~Hhn)  c  Iso...sn-l-  Hence,  the  intersection 
nn>oIso...sn    is  nonempty  for  any  n  £  0.  Note  that  if   x  e  nn^0Iso...sn/  then 
x  €  Iso,  f(x)  €  Isi, . . . ,  so  that  S(x)  =  S0S1S2 This  proves  that  S  is  onto. 

Since  S  is  one-to-one  and  onto,  it  follows  that  S"1  exists.  Hence  we  need 
only  show  that  S  is  continuous  to  prove  it  is  a  homeomorphism.   Thus, 

choose  x  €  A  and  let  S(x)  =  S0S1S2 Let  e  >  0,  and  choose  n  such  that  1  /2n 

<  e.  Consider  It0...tn  for  all  possible  combinations  of  to,  ti, . . . ,  tn.  These  sub- 
intervals  are  disjoint,  and  A  is  contained  in  their  union.  Choose  y  e  A  and 
5  such  that   I  x  -  y  I  <  6.  Then  y  e  Iso...sn/  and  S(x)  and  S(y)  agree  for  the  first 
n+1  terms.  Hence,  by  the  metric  on  code  space  £2/  d(S(x),  S(y))  <  l/2n  <  e. 
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Thus  S  is  continuous.   Trivially,  S"1  is  also  continuous;  hence,  S  is  a 
homeomorphism. 

Since  S:  A— >X  is  a  homeomorphism  ,  we  use  topological  conjugacy  to 
show  that  f  has  the  same  dynamics  as  the  shift  map  a  on  code  space  £2-  The 
commutative  diagram  for  this  relationship  is  shown  in  Figure  4.34.   Since  we 
know  a  is  chaotic  on  X,    f  is  chaotic  on  A  through  topological  conjugacy. 

A      —  f    ->      A 

is  is 

X     —  a  ->     X 

Figure  4.34  The  topological  conjugacy  between  f:  A-»A  and  a:  X-»X. 

We  now  return  to  the  map  f(x)  =  x2  -  1.755  to  learn  more  about  it 
through  its  symbolic  dynamics.  Since  it  has  a  three-cycle,  -1.755,  0,  and 
1.325. . . ,  Sarkovskii's  theorem  guarantees  it  is  chaotic.  However,  we  see  we 
could  have  discovered  this  feature  without  Sarkovskii's  theorem. 

We  start  by  finding  three  open  intervals,  Oi  about  0,  O2  about  -1.755, 
and  O3  about  1.325.  Select  these  such  that  Oi  contains  the  closure  of  f(Oi), 
and  Oi+i  z>  f(Oi).  We  can  always  make  this  choice  because  f  is  a  continuous 
map  (although  in  practice,  the  use  of  a  computer  would  help).  Now  let  Io 
denote  the  closed  interval  between  Oi  and  O3  and  let  Ii  be  the  closed 
interval  between  O2  and  Oj.  This  relationship  is  shown  schematically  in 
Figure  4.35.  We  may  choose  each  Oi  such  that  I  (^'(x)  I  >  1  on  Io  u  Ii.  We 
then  have  f(In)  =>  Ii  and  f(Ii)  z>  Io  U  Ii,  so  each  interval  is  stretched  over  its 
image. 
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-1.755                    0                  1.325 
«•— ( 1 ) ( 1 ) f 1 ) ► 

°2      [    l1    ]       °1       [l2]     °3 

Figure  4.35  The  closed  intervals  I  and  the  open  intervals  O. 


We  now  introduce  symbolic  dynamics.  Let  A  =  {x:  f^x)  e  Io  u  Ii  V  n  >  0). 
We  know  that  A  is  a  Cantor  set.  To  model  the  dynamics  of  f  on  A,  consider 
modified  code  space  X'  where 

X'  =  {S1S2S3. . . :  Si  e  {0, 1}  and  Sk  =  0  =>  Sk+i  =  1}/ 
i.e.,  this  is  just  £2  with  no  adjacent  pairs  of  0s.  If  we  now  define  the  map 
S:  A— >£'  as  above,  we  see  that  the  condition  f(Io)  z>  Ii  forces  the  condition  of 
no  adjacent  0s  in  X*.  The  diagram  in  Figure  4.36  shows  how  f  commutes 
with  a  through  S.  Hence  the  shift  map  a:  £'-»£'  provides  all  of  the 
information  about  the  dynamics  of  f:  A-»A.  Thus  there  exist  points  of  all 
periods  in  £'.  In  fact,  the  point  0111. .  .10111. .  .10. . . ,  with  blocks  of    n-1 
repeating  Is,  is  the  same  point  found  in  the  proof  of  the  special  case  of 
Sarkovskii's  theorem. 


A     —  f    ->     A 

is  is 

I1  —  a  -»     I' 

Figure  4.36  Topological  conjugacy  between  f  and  a. 
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L.    JULIA  SETS 

Julia  sets,  along  with  the  Mandelbrot  set,  have  perhaps  been  the  most 
significant  factors  in  generating  interest  in  chaos  among  laymen.   The  reason 
is  because  the  intricate  and  beautifully  colored  computer  images  shown  as 
"pictures  of  chaos"  are  normally  pictures  of  Julia  sets.   Moreover,  the 
"movies"  of  these  images,  exploding  across  the  screen,  are  simply  the  Julia 
sets  viewed  under  the  continuous  changing  of  a  parameter. 

Julia  sets  were  actually  discovered  in  the  1920s  by  the  French 
mathematicians  Gaston  Julia  and  Pierre  Fatou.   However,  their  true  beauty 
and  intricate  detail  were  not  fully  realized  until  the  1970s  when  computer 
graphics  allowed  for  their  inspection  in  detail.  The  concept  of  a  Julia  set  can 
be  understood  with  only  a  basic  understanding  of  complex  numbers.  On  the 
other  hand,  a  formal  and  mathematically  rigorous  treatment  of  Julia  sets 
requires  a  theory  of  complex  analysis  beyond  the  scope  of  this  thesis.  Here  we 
present  only  a  cursory  survey  of  Julia  sets  in  their  ambient  space,  the  complex 
plane,  still  treating  one-dimensional  maps  in  the  iterated  function  systems. 
While  many  of  the  references  cited  discuss  Julia  sets,  the  presentation  here  is 
from  Keen  (1989)  and  Falconer  (1990).  There  are  many  equivalent  definitions 
of  Julia  sets,  but  the  one  we  present  is  perhaps  the  simplest  to  demonstrate 
and  understand. 

DEFINITION.  Given  a  mapping  f:  C-»C  of  the  complex  plane,  its  Julia  set 
J(f)  is  the  closure  of  the  set  of  repelling  periodic  points  of  f. 


131 


As  a  simple  example,  consider  the  map  f(z)  =  2z.  Under  infinite  iteration 
of  f,  all  points  in  the  complex  plane  excluding  the  origin  tend  to  *»  (or  to  be 
more  precise,  the  point  at  <*>).  The  origin  is  a  fixed,  and  hence  periodic  point 
of  this  map.  Since  iteration  of  any  point  other  than  the  origin  tends  away 
from  it,  the  origin  is  repelling.  This  result  is  also  verified  from    I  f  (0)  I  =  121  > 
1.  Since  the  only  repelling  periodic  point  of  this  map  is  the  origin,  which  is 
its  own  closure,  the  Julia  set  for  f(z)  =  2z  is  the  origin. 

We  now  present  a  less  trivial  example  which  demonstrates  many 
interesting  properties  of  Julia  sets.  Consider  the  map  f (z)  ■  z2  +  c  for  c  ■  0. 
All  points  inside  the  unit  circle     I  z  I  <  1  tend  to  the  origin  under  infinite 
iteration.  Thus  the  origin  is  an  attracting  fixed  point  of  the  map.  In  fact, 
I  f  (0)  I  =  1 2(0)  I  =  0,  which  also  verifies  that  the  origin  is  attracting.  Moreover, 
all  points    I  z  I  >  1  outside  the  unit  circle  tend  to  ~  under  iteration  of  this 
map. 

Now  consider  the  standard  unit  circle,  I  z  I  =1.  These  points  are 
represented  by  z  =  eie.  Then  z2  =  e12®,  which  is  exactly  the  chaotic  map  of  the 
unit  circle  f:  S1-»S1,  where  f(0)  =  20,  studied  earlier.  We  know  the  periodic 
points  of  this  map  are  dense  on  the  unit  circle.  Since  the  periodic  points  of 
f  (z)  =  z2  are  dense  on  the  unit  circle,  every  point  on  the  unit  circle  is  the  limit 
of  a  sequence  of  periodic  points  of  f .  Thus,  the  closure  of  the  periodic  points 
of  f  is  the  unit  circle.  Moreover  these  points  are  repelling.  To  see  this,  recall 
that  points  inside  the  unit  circle  converge  to  the  origin,  whereas  points 
outside  the  unit  circle  diverge  to  the  point  at  infinity.  Furthermore, 
I  f  (z)  I  1 1  z  |  -i  =  1 2(1)  I  =  2  >  1,  verifying  that  these  points  are  repelling.  Hence, 
the  Julia  set  for  the  function  f(z)  =  z2  is  the  standard  unit  circle. 
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Another  concept  associated  with  a  Julia  set  is  that  of  the  filled  Julia  set, 
denoted  F(f).  When  the  Julia  set  is  a  closed  curve,  the  set  F(f)  is  the  union  of 
J(f)  with  its  interior.  The  filled  Julia  set  is  the  set  of  points  that  do  not  escape 
to  oo  under  infinite  iteration  of    f.  For  the  example  f(z)  =  z2,  the  filled  Julia 
set  is  the  closed  unit  disc,  F(f)  =  {z:    I  z  I  <  1).  As  an  aside,  the  complement  of 
the  Julia  set  is  called  the  Fatou  set  and  is  sometimes  denoted  F(f)  as  well, 
although  Jc  is  also  used.  Loosely  speaking,  J(f)  is  the  set  containing  the 
"bad"  (i.e.,  chaotic)  behavior,  while  the  Fatou  set  is  the  "good"  set,  possesing 
the  well-behaved  dynamics. 

Having  introduced  the  concepts  of  Julia  sets  and  filled  Julia  sets  in  this 
simple  setting,  we  now  describe  an  algorithm  for  generating  computer  images 
of  these  objects.  If  we  superimpose  the  complex  axes  on  a  computer  screen  to 
an  appropriate  scale,  then  points  in  the  complex  plane  correspond  to  pixels  on 
the  screen,  although  this  relationship  is  certainly  not  one-to-one.  Given  a 
function  f(z),  we  can  iterate  each  pixel.  Since  we  are  interested  in  the  points 
that  escape  to  <*>,  a  bound  (normally  very  large)  can  be  set  which  we  call   I Z I , 
and  above  which  an  iterate  is  considered  as  having  escaped.  Next  select  k 
integers    N\  <  N2  <  ...  <  N^  =  N.  Color  the  screen  with  k  +  1  colors  based 
on  the  following  algorithm:  as  a  point  is  iterated,  if  it  has  not  escaped  after  N 
iterations,  color  it  black.  If  it  escapes  (goes  beyond  I Z I )  between  0  and  N\ 
iterations,  assign  to  it  another  color  (say,  red).  If  it  escapes  between  Ni  and 
N2  iterations,  color  it  with  yet  another  color  (for  example,  yellow).  Continue 
in  this  manner  until  the  entire  screen  has  been  colored.  Selecting  a  large 
value  for  k  provides  more  detail,  which  can  be  refined  further  by 
experimentally  adjusting  N  and    I Z I   with  respect  to  each  other.  The  part  of 
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the  screen  colored  black  (if  we  have  chosen  N,    I Z I ,  and  the  scale 
appropriately)  is  the  filled  Julia  set  for  f,  and  this  region's  boundary  (provided 
it  is  connected)  is  the  Julia  set  itself.  The  colored  bands  around  the  Julia  set 
are  contours  corresponding  to  the  various  escape  times  of  points  in  the 
exterior  of  the  Julia  set.  The  reason  the  complement  of  the  filled  Julia  set  is 
colored  is  because  of  the  finite  scale  of  the  computer  screen:  there  can  be  great 
detail  occurring  within  the  area  of  a  single  pixel  and,  while  the  complete  Julia 
set  is  not  revealed  by  just  the  black  area,  much  can  be  determined  about  its 
border  by  examining  the  distorted  contours  surrounding  it. 

For  the  example  f(z)  =  z2,  coloring  the  screen  based  on  this  escape  time 
algorithm  produces  a  black  disc  with  a  sequence  of  concentric  colored  circles 
around  it  (which  in  itself  is  not  particularly  interesting).  However,  recalling 
the  family  of  functions  fc(z)  =  z2  +  c,  as  the  parameter  c  is  varied  some  very 
interesting  results  occur.  Unfortunately,  while  this  Julia  set  has  many 
fascinating  properties,  an  advanced  level  of  complex  analysis  is  required  to 
establish  even  its  most  basic  properties.  The  required  concepts  include 
families  of  normal  functions,  the  Arzela-Ascoli  theorem,  and  Montel's 
theorem.  These  results  are  beyond  the  scope  of  this  thesis,  but  an  excellent 
summary  of  them  is  found  in  Falconer  (1990),  and  we  provide  a  synopsis  of 
them  at  the  end  of  this  section.  Nevertheless,  we  can  still  provide  a  brief 
description  of  some  of  the  salient  characteristics. 

As  the  parameter  c  is  varied  away  from  the  origin,  the  Julia  set  (the  unit 
circle)  begins  to  continuously  distort  and  take  on  different  shapes.  Closer 
inspection  reveals  that  the  boundary  appears  to  become  infinitely  detailed 
and  self-similar;  in  fact,  it  becomes  fractal.  (Even  with  c  =  0,  the  boundary  of 
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the  filled  Julia  set  is  infinitely  self-similar,  although  in  a  trivial  manner.) 
While  we  saw  that  the  dynamics  of  f  on  J(f)  for  c  =  0  are  chaotic  (since  they 
share  the  properties  of  f(6)  =  29)  they  continue  to  be  chaotic  on  J(f)  as  c 
varies.  For  the  quadratic  family  f(z)  =  z2  +  c,  as  c  is  varied  the  Julia  set  varies 
from  a  circle,  to  a  closed  curve  with  "bulbs"  that  are  "pinched"  together  at  a 
single  point,  to  "dendrites"  which  are  fractal  structures  with  no  interior,  to 
"dust"  which  is  a  set  of  disconnected  points  which  are  scattered  about  a  region 
of  the  complex  plane,  similar  to  the  Cantor  set.  For  some  Julia  sets  with 
fractal  boundaries,  like  the  Koch  snowflake,  the  lengths  of  the  boundaries  are 
infinite.  A  further  result  about  Julia  sets  is  that  they  are  either  connected 
(meaning  they  consist  of  one  solid  piece)  or  totally  disconnected  (meaning 
they  have  a  structure  similar  to  Cantor  dust).  The  Julia  sets  for  various 
values  of  c  are  shown  in  Figure  4.37. 

Without  going  into  too  much  detail,  we  provide  a  brief  synopsis  of  the 
most  important  ingredients  of  the  mathematical  theory  behind  these  results. 
For  analytic  functions  in  C  (i.e.,  those  that  are  infinitely  differentiable  in  the 
complex  sense)  techniques  of  complex  variable  theory  can  be  used  to  establish 
the  basic  properties  of  Julia  sets. 

It  can  be  shown  that  an  alternative  (but  equivalent)  way  of  defining  the 
Julia  set  J(f)  for  polynomials  f  is  as  the  set  of  all  complex  values  z  for  which 
the  family  {fk(z)}  k  =  1,  2, . . .  is  not  normal.  Loosely  speaking,  a  family  of 
complex  analytic  functions  is  said  to  be  normal  if  it  possesses  some  expecially 
strong  convergence  properties  (called  "uniform"  convergence)  on  compact 
subsets  of  a  given  open  set.  Using  a  powerful  result  from  complex  analysis 
known  as  Montel's  theorem,  it  is  possible  to  demonstrate  that  if  f  is  a 
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Figure  4.37  Julia  sets  of  the  function  f(z)  =  z2  +  c  for  (a)  c  =  -.1  +  .li; 

(b)  c  =  -.5  +  .5i;  (c)  c  =  -1  +  .05i;  (d)  c  =  -.2  +  .751;  (e)  c  =  .25  +  .52i;  (f)  c  = 

-.5  +  .55i;  (g)  c  =  .66i;  (h)  c  =  -i.  The  figure  is  from  Falconer  (1990,  p.  213). 
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polynomial,  then  J(f)  is  non-empty,  compact  (closed  and  bounded),  contains 
no  "isolated"  points,  and  much  more.   Note  that  this  is  consistent  with  what 
we  have  seen  with  the  quadratic  family,  even  in  the  case  when  J(f)   is  "dust." 
It  should  be  pointed  out  that  the  above  results  do  not  necessarily  hold  for 
non-polynomial  complex  analytic  maps.   The  Julia  set  for  the  exponential 
map  f(z)  =  e2,  for  example,  is  the  entire  complex  plane.  Of  course, 
polynomial  functions  are  not  the  only  ones  that  generate  interesting  Julia 
sets.  Some  of  the  trigonometric  families,  such  as  X  sin  z,  also  provide  very 
interesting  characteristics  as  X,  varies. 

M.  THE  MANDELBROT  SET 

The  Mandelbrot  set  is  often  associated  with  intricate  computer  graphics. 
It  has  been  described  from  "the  most  complex  object  in  mathematics"  to  "the 
most  beautiful  object  in  mathematics."  While  Julia  sets  are  found  in  range 
space  of  a  complex  function,  the  Mandelbrot  set  lies  in  parameter  space, 
which  is  the  complex  plane  when  the  parameter  is  a  complex  number.  Like 
the  Julia  set,  almost  every  reference  cited  discusses  the  Mandelbrot  set. 
However,  th  particular  presentation  here  is  based  on  Branner  (1989)  and 
Falconer  (1990).  There  are  two  equivalent  definitions  of  the  Mandelbrot  set, 
and  both  of  them  are  presented  here. 

One  definition  of  the  Mandelbrot  set  for    fc(z)  =  z2  +  c  is  the  set  of  values 
of  c  for  which  the  associated  Julia  set  J(f)  is  connected.  (This  definition 
stresses  the  connection  between  the  Mandelbrot  set  and  Julia  sets.)  As 
mentioned  in  the  last  section,  the  Julia  sets  for  f(z)  =  z2  +  c  vary  from  being 
totally  connected  to  "dust."  The  values  of  c  for  which  the  Julia  sets  are  dust 
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do  not  belong  to  the  Mandelbrot  set.   An  equivalent  definition,  which  is 
perhaps  easier  to  understand,  follows. 

DEFINITION.  The  Mandelbrot  set  Tt  is  the  set  of  complex  values  of  c  for 
which  the  origin  does  not  escape  to  ©°  under  infinite  iteration  of  f(z)  =  z2  +  c. 

A  picture  of  the  Mandelbrot  set  in  the  complex  plane  is  shown  in  Figure 
4.38.  The  figure  is  from  Falconer  (1990,  p.  205).  We  have  already  seen  that  for 
c  >  1/4  on  the  real  line,  all  values  of  x  including  the  origin  go  to  infinity 
under  iteration  of  f(x)  =  x2  +  c;  for  c  <  -2,  the  origin  also  escapes.  Hence,  we 
know  that  the  Mandelbrot  set  contains  the  interval  [-2, 1/4]  on  the  real  line. 
However,  the  situation  is  not  nearly  this  simple  when  c  varies  in  the 
complex  plane. 
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Figure  4.38  The  Mandelbrot  set. 
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Like  Julia  sets,  the  Mandelbrot  set  is  infinitely  detailed.  In  fact,  it  contains 
many  smaller  copies  of  itself  around  its  border.  However,  it  is  by  no  means 
self-similar  because  it  contains  many  other  interesting  shapes.   Moreover,  its 
intricate  detail  varies  significantly  among  the  border  regions  of  the  smaller, 
Mandelbrot-like  sets. 

Pictures  of  the  Mandelbrot  set  can  also  be  generated  using  the  escape  time 
algorithm  used  to  draw  Julia  sets.  Here  the  coloring  of  the  complement 
becomes  particularly  important  because  many  of  the  tendrils  extending  from 
the  main  body  of  the  Mandelbrot  set  are  too  detailed  to  capture  on  a  computer 
screen  (regardless  of  the  scale  chosen),  so  they  are  only  evidenced  by  the 
distorted  contours  surrounding  them.  It  is  known  that  the  Mandelbrot  set  is 
connected:  even  points  that  appear  isolated  on  computer  images  are 
connected  to  the  main  body  by  dendrites  too  small  to  be  seen  on  a  computer 
screen. 

The  first  definition  of  the  Mandelbrot  set  shows  an  intimate  connection 
with  Julia  sets  of  the  function  f(z)  =  z2  +  c,  but  the  Mandelbrot  set  contains 
even  more  information  about  the  dynamics  of  the  function  f .   Once  again, 
however,  any  rigorous  mathematical  development  of  these  dynamics 
requires  the  advanced  theory  of  complex  analysis.  So  again  we  only  describe 
some  of  the  more  interesting  characteristics. 

The  first  result  is  that  each  "bulb"  of  the  Mandelbrot  set  corresponds  to  an 
attracting  k-cycle  of  f (z)  for  a  particular  value  of  k.  For  example,  the  large 
central  cardioid  corresponds  to  the  values  of  c  for  which  f(z)  =  z2  +  c  has  an 
attracting  fixed  point.  To  see  this,  note  that  an  attracting  fixed  point  must 
satisfy   z2  +  c  =  z  and   lf(z)l  =  I2zl  <  1.  The  boundary  of  this  region  is  given 
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by  c  =  z  -  z2,  where  I z I  =1/2.  In  polar  representation,  this  becomes 

c  =  (1/2)  e2™9  -  (1/4)  e4*19,  0  <  0  <  2k.  These  values  of  c  trace  out  a  cardioid 

in  the  complex  plane  with  a  cusp  at  z  =  1/4  +  Oi.  Unfortunately,  the  periods 

of  the  attractive  cycles  of  the  other  bulbs  do  not  so  easily  reveal  themselves 

mathematically. 

It  is  perhaps  not  surprising  that  the  periods  of  the  bulbs  along  the  real  axis 
are  in  direct  correspondence  with  the  bifurcations  found  for  the  map 
f (x)  =  x2  +  c.  Recall  that  for  this  latter  map,  a  tangent  bifurcation  occurs  at 
c  =  1/4,  and  a  series  of  period-doubling  bifurcations  begins  as  c  decreases 
through  -3/4.  Figure  4.39  shows  the  Mandelbrot  set  plotted  on  the  same 
coordinate  axis  as  this  bifurcation  diagram.  You  can  see  the  alignment  of  the 
main  bulbs  with  the  period  doubling  that  occurs  along  the  real  axis.  The  bulb 
in  the  "tail"  of  the  Mandelbrot  set  corresponds  to  the  three-cycle  that  emerged 
out  of  chaos  around  the  value  c  ■  -1.755  studied  earlier. 

The  Julia  sets  associated  with  the  c  values  belonging  to  the  Mandelbrot 
set  vary  as  the  period  of  the  attracting  cycle  varies  among  the  bulbs.  Julia  sets 
for  values  of  c  in  some  of  the  different  bulbs  of  the  Mandelbrot  set  are  shown 
in  Figure  4.40.  Notice  that  the  number  of  "bulbs"  in  the  Julia  sets  that  are 
pinched  together  at  a  single  point  correspond  to  the  period  of  the  cycles  of  the 
Mandelbrot  set.  For  example,  the  Julia  sets  for  values  of  c  in  the  main 
cardioid  are  all  simple  closed  curves  which  correspond  to  the  attractive  fixed 
points,  whereas  the  values  of  c  in  bulbs  that  correspond  to  attractive  n-cycles 
have  n  bulbs  converging  at  a  single  point.  Notice  also  the  very  thin  Julia  set 
(dendrite)  associated  with  one  of  the  tendrils  of  Tt.  Dendrites  occur  for 
values  of  c  for  which  the  origin  is  a  periodic  point  of  f (z)  =  z2  +  c;  for 
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Figure  4.39  The  Mandelbrot  set  plotted  against  the  bifurcation  diagram  for 
f(x)  =  x2  +  c.  The  figure  is  from  the  article  by  Devaney  (1989,  p.  37). 


example  the  point  c  =  -i  as  shown  in  Figure  4.37.  Finally,  the  Julia  set 
associated  with  a  point  not  in  the  Mandelbrot  set  is  totally  disconnected. 
Remembering  that  for  values  of    c  <  -2,  the  sets  of  periodic  points  for  the 
iterated  maps  f(x)  =  x2  +  c  were  Cantor  sets  (hence  totally  disconnected)  their 
Julia  sets,  by  definition,  are  also  disconnected.  Thus  these  points  fail  to  belong 
to  the  Mandelbrot  set. 


141 


*& 


M 


~k 


^M^ 

♦  ••V**  5 


Figure  4.40  The  Julia  sets  of  points  in  different  bulbs  of  the  Mandelbrot  set 
(see  Figure  4.37).  The  figure  is  from  Falconer  (1990,  p.  214). 


Another  feature  of  the  Mandelbrot  set  is  the  existence  of  a  dense  set  on  its 
boundary  of  points,  called  Misiurewicz  points,  for  which  the  image  of  the 
Mandelbrot  set  in  parameter  space,  and  the  corresponding  Julia  set  in  the 
range  space,  look  the  same  up  to  a  rotation  (in  a  sense  that  can  be  made 
mathematically  precise;  see  Branner  (1989)).  Figure  4.41  shows  a  blowup  of 
the  Mandelbrot  set  and  the  Julia  set  around  the  Misiurewicz  point 
c«  -.101096  +  K.956287). 

The  Mandelbrot  set  occurs  in  spaces  other  than  the  parameter  space  we 
have  presented.  In  fact,  it  appears  to  be  an  almost  universal  geometric  shape. 
Recall  the  coloring  of  the  complex  plane  through  Newton's  method  for  the 
function  z4  -  1  =  0  in  Section  IVG.  If  we  color  the  complex  plane  for  different 
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Figure  4.41   The  Julia  set  and  the  Mandelbrot  set  around  a  Misiurewicz  point. 
The  figure  is  from  Branner  (1989,  p.  103). 

values  of  the  parameter  X  for  cubic  polynomials  of  the  family  Px,  chaotic 
regions  are  found  between  the  basins  of  attraction.  However,  interspersed 
within  these  regions  are  small  copies  of  the  Mandelbrot  set.  While  these 
regions  have  always  existed,  it  has  taken  present-day  powerful  computer 
graphics  to  reveal  them.  As  scientists  continue  to  use  computers  to  examine 
dynamical  systems  more  closely,  we  expect  that  the  Mandelbrot  set  will 
appear  with  increased  frequency. 

There  are,  no  doubt,  other  more  fascinating  properties  of  the  Mandelbrot 
set  yet  to  be  discovered  or  proven.  Each  property  reveals  something  about  the 
complexity  of  the  iterated  map  f(z)  =  z2  +  c.  While  we  have  not  discussed  all 
known  results  here,  this  cursory  summary  does  provide  considerable  insight 
into  the  complexity  of  this  chaotic  mapping. 

While  the  function  f(z)  =  z2  +  c  appears  to  be  a  very  specific  form  of  the 
quadratic  family,  it  is  in  fact  topologically  conjugate  to  every  quadratic 
function  for  various  values  of  c.  To  see  this,  consider  the  function 
H(z)  =  az  +  p  with  a  *  0.  Then   h-l(f(h(z))  =  (a2z2  +  2apz  +  p2  +  c  -  p)/a. 
Appropriate  choices  for  the  values    a,  p,  and  c  produce  any  quadratic 
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function  whatsoever.  Thus,  in  studying  the  dynamics  of    f (z)  -  z2  +  c,  reveals 
the  entire  family  of  quadratic  functions. 

N.  THE  SMALE  HORSESHOE 

Our  attention  so  far  has  been  restricted  to  one-dimensional  dynamical 
systems.  In  so  doing,  we  have  learned  much  about  chaos.  However,  since 
many  real-world  phenomena  occur  in  two  and  three  dimensions,  the  range 
of  applications  has  been  restricted.  We  now  investigate  our  first  two- 
dimensional  dynamical  system,  the  Smale  horseshoe.   Instead  of  developing 
the  horseshoe  algebraically,  a  strict  geometric  interpretation  of  the  map  is 
given.  The  primary  references  for  this  section  are  Holmes  (1989), 
Guckenheimer  (1990),  and  Devaney  (1989). 

The  Smale  horseshoe  was  originally  constructed  to  help  interpret  the 
periodically  forced  oscillator,  which  commonly  appears  in  applications  in 
physics,  mechanics,  and  electrical  engineering.  Normally,  the  systems  under 
investigation  are  modeled  with  ordinary  differential  equations,  and  the 
Smale  horseshoe  turns  out  to  provide  an  intuitive  way  to  see  why  the 
equations  sometimes  lead  to  chaotic  behavior. 

Many  versions  of  the  Smale  horseshoe  exist.  We  present  here  the 
version  that  is  the  simplest  geometrically.  Thus,  take  the  unit  square  in 
Figure  4.42,  stretch  it  out  by  a  factor  of  three  in  one  direction,  and 
simultaneously  shrink  it  by  a  factor  of  three  in  the  other  direction  to  obtain  a 
long  bar.  Then  bend  the  middle  section  of  the  bar  into  a  horseshoe  and 
superimpose  it  back  on  the  original  square,  as  shown  in  the  figure.  Denote 
this  geometric  mapping  by  F.  Notice  that  the  two  shaded  bands  do  not  escape 
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the  unit  square  under  this  first  iteration.   Their  preimages  are  the  horizontal 
bands  shown  in  the  figure. 

Because  the  preimage  of  F  can  be  determined  precisely,  it  is  an  invertible 
map.  Thus  it  is  possible  to  study  not  only  the  forward  orbit  of  points,  but 
their  backward  orbits  as  well.  We  are  interested  in  finding  the  invariant  set  of 
the  unit  square  under  the  forward  and  backward  orbits  of  F.  These  are  the 
points  which  do  not  escape  the  unit  square  under  infinite  forward  and 
backward  iteration.  Then  we  will  be  able  to  investigate  the  dynamics  of  the 
particular  physical  system  associated  with  the  Smale  horseshoe  by  studying 
the  dynamics  on  this  invariant  set. 
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Figure  4.42  Construction  of  the  Smale  horseshoe. 
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Now,  iterate  the  map  a  second  time,  as  shown  in  Figure  4.43.  Observe 
that  the  image  of  the  shaded  area  of  the  first  iteration,  and  its  preimage  appear 
as  before. 
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Figure  4.43  The  second  iteration  of  the  Smale  horseshoe. 


By  superimposing  the  image  of  F  on  its  preimage  for  the  first  two  iterates,  we 
construct  geometrically  an  invariant  set,  shown  as  the  darkly  shaded  region 
in  Figure  4.44.  Here  we  label  the  horizontal  and  vertical  bands  H  and  V, 
respectively. 
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Figure  4.44  The  invariant  set  of  the  Smale  horseshoe. 

Notice  that  the  set  of  points  A  remaining  in  the  unit  square  under 
infinite  forward  and  backward  iteration  has  a  Cantor-like  appearance.  In  fact, 
that  set  turns  out  to  be  the  direct  product  of  two  Cantor  middle-thirds  sets. 
The  variations  of  the  Smale  horseshoe  mentioned  earlier  involve  using 
different  values  for  shrinking  and  stretching  the  unit  square  under  F,  as  well 
as  using  a  different  placement  of  the  horseshoe  when  it  is  superimposed  back 
on  the  square.   All  variations,  however,  still  create  Cantor-like  invariant  sets. 

In  order  to  understand  the  dynamics  of  this  system,  we  only  need  to 
analyze  the  dynamics  on  the  invariant  set  (since  all  other  points  escape  under 
iteration  for  F).  To  undergo  this  analysis,  first  note  that  the  forward  and 
backward  orbits  of  any  point  x  in  the  invariant  set  also  belong  to  it. 
Specifically,  each  point  in  these  orbits  is  in  one  of  the  horizontal  bands  Ho  or 
Hi.  Hence,  define  the  mapping  S:  A-»X  by  the  rule  Sj(x)  =  i  if  Fj(x)  €  Hi  for 
i  e  {0, 1}.  Thus,  every  point  x  in  the  invariant  set  is  associated  with  an 
infinite  string  of  indices  of  the  horizontal  bands  to  which  it  is  mapped  under 

F.  Notice  that  the  index  of  points  in  X  runs  j  =  . . . ,  -3,  -2,  -1,  0, 1, 2, 3, So 

unlike  code  space  X2  for  one-dimensional  maps  (which  consisted  of  semi- 
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infinite  sequences),  the  space  X  consists  of  bi-infinite  sequences.  Figure  4.45 
shows  a  point  x  and  its  orbit  under  three  forward  iterates  and  one  backward 
iterate.  So  for  j  =  . . . ,  -1,  0,1,2,3,...,  Sj(x)  =  ...  00100. . .  from  the 
horizontal  bands  in  which  each  iterate  lies.  Noticing  that  Sj(F(x))  =  Sj+i(x), 
one  sees  that  F  applied  to  the  set  A  corresponds  to  the  shift  map  a  on  the 
space  £.  Moreover,  every  symbol  in  £  corresponds  to  a  unique  orbit  of  F, 
because  every  image  V  completely  intersects  its  preimage  H.  Therefore,  the 
mapping  F  and  the  shift  map  o  on  infinite  code  space  are  topologically 
conjugate  through  the  map  S,  as  shown  in  Figure  4.46. 

The  horseshoe  map  has  been  very  useful  in  analyzing  physical  systems 
because  it  extends  to  any  Euclidean  space  Rn.  The  connection  with  ordinary 
differential  equations  is  through  a  concept  known  as  the  Poincare  map.  If  the 
phase  space  associated  with  an  ordinary  differential  equation  is  intersected 
with  a  plane  normal  to  any  orbit,  then  the  orbit  intersects  the  plane  exactly 
once  during  each  cycle.  The  collection  of  these  points  of  intersection  is  called 
the  Poincare  map.  While  the  horseshoe  map  was  constructed  originally  in 
connection  with  the  Poincare  map  of  a  periodically  forced  oscillator,  there  is  a 
general  method  for  finding  horseshoes  that  applies  to  a  wide  range  of 
Poincare  maps.  The  procedure  has  helped  scientists  and  engineers 
understand  the  dynamics  of  the  associated  physical  systems. 

The  actions  of  stretching  and  bending  in  the  Smale  horseshoe  are 
frequently  encountered  in  physical  systems.  Predicting  the  orbit  of  points  in 
such  systems  (a  simple  taffy  pull  serves  as  a  classical  example)  has  always 
proven  elusive.  The  science  of  chaos  has  helped  explain  why  these  systems 
have  been  so  difficult  to  understand. 
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Figure  4.45  The  orbit  of  a  point  x  of  the  invariant  set 
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X     —  o  -»    I 
Figure  4.46  The  topological  conjugacy  between  F:  A-»A  and  a:  X-»£. 


O.    THEHENONMAP 

With  the  Smale  horseshoe  providing  a  geometric  introduction  to  two- 
dimensional  dynamical  systems,  we  now  turn  our  attention  to  another  map 
of  the  plane  that  exhibits  many  of  the  interesting  properties  of  two- 
dimensional  maps.  The  material  in  this  section  is  presented  as  a  series  of 
exercises  in  Devaney  (1989),  to  which  most  of  the  answers  and  results  come 
from  Rasband  (1990),  Alligood  (1989),  Cherbit  (1991)  and  Moon  (1987). 
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The  Henon  map  H^:  R2-»R2,  is  defined  by  the  equations 

xi  =  1  +  yo  -  axo2, 

yi  =  bxo. 

Notice  that  H  depends  on  two  parameters,  a  and  b,  and  that  it  has  only  one 
nonlinear  term   (x2).  Thus  H  is  one  of  the  simplest  higher-dimensional 
nonlinear  maps  we  can  study.  A  number  of  questions  regarding  the  Henon 
map  have  not  been  resolved  because  of  the  wide  range  of  possible  parameter 
values,  but  for  certain  parameter  values  it  exhibits  some  very  interesting 
behavior. 

First,  note  that  the  Henon  map  can  be  expressed  as  the  composition  of 
three  maps  H3°H2°Hi,  where  Hi(x,  y)  =  (x,  1  -  ax2  +  y)  is  a  nonlinear  bending 
(and  a  quick  check  with  calculus  shows  it  is  area  preserving);  H2(x,  y)  =  (bx,  y) 
is  an  expansion  or  contraction  in  the  x  direction,  depending  on  the  value  of 
b;  and  H3(x,  y)  ■  (y,  x)  flips  the  contracted,  bent  image  about  the  main 
diagonal. 

The  case  b  =  0  makes  the  Henon  map  topologically  conjugate  to  the  map 
g(x)  =  1  -  ax2  if  we  consider  the  projection  of  H  onto  the  x-axis.  For  the  case 
I  b  I  >  1,  the  map  H2  is  not  a  contraction  and  the  iterates  diverge.  Hence,  we 
restrict  our  attention  to  the  range  0  <  I  b  I  <  1. 

Now  fix  b.  It  is  easy  to  show  the  fixed  points  of  the  Henon  map  are 
(x,  y)  =  ([b  - 1  ±  ([b  - l]2  +  4a)V2]/2a,  bx). 
A  closer  inspection  reveals  that  for  (b  - 1)2  +  4a  <  0,  or  a  <  -(b  - 1)2/4,  these 
points  have  an  imaginary  component  yielding  no  fixed  points  in  R2. 
Moreover,  when  a  =  -(b  -  l)2/4,  the  fixed  points  coincide  (so  there  is  only  one 
attracting  fixed  point).  Finally,  for  a  >  -(b  -  l)2/4,  there  are  two  distinct  fixed 
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points,  one  of  which  is  attracting.  Here,  for  a  fixed  value  of  b,  as  the 
parameter  a  increases  through  a  critical  parameter  value,  a  tangent 
bifurcation  occurs. 

As  the  parameter  a  continues  to  increase  a  series  of  period-doubling 
bifurcations  appears  eventually  leading  to  chaos.  Let  a*,  denote  the  value  of 
a  beyond  which  chaos  occurs.  Then  the  dynamics  of  the  Henon  map  can  be 
determined  geometrically  in  a  familiar  setting.  For  a  fixed  value  of  b,  let  R 
be  the  larger  root  of  a£2  -  (b  -  1)£  -  1  -  0.  Let  S  be  the  square  centered  at  the 
origin  with  vertices  at  (±R,  ±R).  Figure  4.47  shows  the  images  of  S  under  H 
for  a  <  aoo  and  a  >  a,*.  Note  also  the  effects  of  Hi,  H2,  and  H3  in  the  way  the 
square  S  is  bent,  contracted,  and  flipped.  Additionally,  for  a  >  a^,  the 
geometric  construction  looks  similar  to  the  Smale  horseshoe  (and,  in  fact,  it  is 
a  horseshoe).  Thus  the  dynamics  of  the  map  H  for  a  >  a<»  (for  a  fixed  b)  are 
indeed  chaotic. 


►  x 


Figure  4.47  The  images  of  S  under  the  Henon  map. 
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The  particular  value  b  =  0.3  has  been  studied  extensively  and  the  tangent 
bifurcation  occurs  at  a  =  .1225.  If  a  increases  holding  b  =  0.3  constant,  in  the 
range  1.052  <  a  <  1.082  a  series  of  period-doubling  bifurcations  occurs  which 
eventually  lead  to  chaotic  behavior. 

A  particularly  interesting  phenomenon  occurs  close  to  b  =  0.3  and    a  -  1.4. 
Here  we  have  an  attractor  of  the  system.  The  infinite  iteration  of  bending, 
shrinking,  and  flipping  the  plane  yields  results  not  yet  fully  understood. 
Nevertheless,  with  the  aid  of  computers,  it  has  been  possible  to  compute 
these  results  numerically  and  view  them  graphically.  Iterating  an  initial 
point  (xo,  yo)  under  H  yields  a  set  of  points,  called  the  attractor  of  H,  that 
appear  to  be  invariant  under  infinite  iteration  of  H.  The  attractor  of  H  for 
the  values  b  =  .3  and  a  ■  1.4  is  shown  in  Figure  4.48.  The  dynamics  on  this 
attractor  are  chaotic  (as  just  shown  geometrically  with  the  analog  to  the  Smale 
horseshoe).  Numerically  it  has  been  found  that  the  attractor  appears  to  have 
a  dense  orbit,  sensitive  dependence  on  initial  conditions,  and  to  be 
topologically  transitive.   However,  since  the  evidence  of  this  invariant  set  has 
only  been  suggested  by  the  use  of  numerical  computation  (and  not  established 
with  any  mathematical  rigor)  many  of  its  properties  are  still  not  clearly 
identified. 

The  Henon  attractor  (if  it  truly  exists)  for  b  =  .3  and  a  =  1.4  fits  into  a 
class  of  attractors  referred  to  as  strange  attractors.  While  a  formal  definition 
of  strange  attractor  has  notbeen  developed  to  date,  there  are  three  conditions 
that  seem  to  be  characteristic  of  them.  These  characteristics  are: 
1.  Points  "nearby"  the  attractor  converge  to  the  attractor  under  infinite 

iteration  of  the  function. 
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2.  The  dynamics  of  points  on  the  attractor  are  chaotic. 

3.  The  attractor  has  a  non-integer  fractal  dimension. 


)      \ 


Figure  4.48  The  Henon  attractor.  The  figure  is  from  Holden  (1986,  p.  90). 

By  "nearby"  we  refer  to  a  region  of  the  plane  called  the  basin  of  attraction, 
inside  of  which  all  points  converge  to  the  attractor.  The  basin  of  attraction 
depends  on  the  particular  function,  but  in  the  case  of  the  Henon  attractor  it 
turns  out  to  be  the  entire  Euclidean  plane. 

Magnification  of  the  Henon  attractor  indicates  that  it  is  infinitely  detailed, 
as  evidenced  by  the  "bands"  in  Figure  4.49  actually  being  composed  of  smaller 
bands  of  points.  Additionally,  its  fractal  dimension  has  been  estimated 
numerically  at  1.26  for  the  parameter  values  b  =  .3  and    a  =  1.4. 
Nevertheless,  considerable  mystery  remains  concerning  the  Henon  attractor 
(as  well  as  many  of  the  other  interesting  strange  attractors  that  have  been 
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discovered  numerically  or  physically).   Because  of  their  structure  and  self- 
similarity,  fractal  geometry  is  currently  being  applied  to  the  study  of  strange 
attractors. 
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Figure  4.49  Magnification  of  the  Henon  attractor.  The  figure  is  from 
Berge(1984,p.  133). 


While  strange  attractors  come  up  in  models  of  physical  equations  such  as 
the  Duffing  equation,  the  van  der  Pol  equations,  or  the  Rossler  equations, 
they  have  also  been  seen  in  physical  systems.  While  many  infectious  diseases 
appear  to  follow  definite  cycles,  measles  appears  to  follow  a  strange  attractor 
with  fractal  dimension  2.5  when  viewed  in  the  proper  phase  space. 
Additionally,  Saturn's  rings,  because  of  their  remarkable  resemblance  to  the 
strange  attractors  of  many  mathematical  systems,  are  being  studied  in  this 
new  light  in  great  detail  (however,  this  connection  is  still  being  investigated, 
and  no  conclusions  have  yet  been  drawn). 

A  final  remark  about  the  Henon  map:  if  we  set  the  parameter  b  =  1,  the 
map  becomes  an  area  preserving  map  of  the  plane.  Since  the  map  Hi(x,y)  = 
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(x,  1  -  ax2  +  y)  is  area  preserving  as  observed  earlier,  we  see  that  with  b  =  1, 
H2(x,  y)  =  (x,  y)  and  H3(x,  y)  =  (y,  x)  also  preserve  areas.  This  condition  leads 
to  an  entirely  new  set  of  phenomena,  one  of  which  we  briefly  mention  here. 
As  the  parameter  a  increases,  orbits  of  different  periods  are  created,  but  the 
last  orbit  to  develop  is  a  two-cycle.  This  provides  an  example  of  where 
Sarkovskii's  theorem  fails  to  apply  in  two  dimensions. 

P.    THE  LORENZ  EQUATIONS 

It  is  appropriate  to  conclude  our  mathematical  treatment  of  chaos  with 
the  Lorenz  equations  because  they  comprise  one  of  the  first  systems  to  bring 
chaotic  dynamical  systems  to  the  attention  of  the  mathematical  community. 
The  primary  references  for  this  section  are  Sparrow  (1982),  Holden  (1986),  and 
Fischer  (1985),  although  some  of  the  presentation  follows  that  of  Berge  (1984), 
Guckenheimer  (1990),  and  Thompson  (1989). 

The  Lorenz  equations  have  been  studied  extensively  since  the  mid  1970s, 
and  numerous  interesting  results  have  been  derived  from  them.   However, 
to  discuss  many  of  these  results  requires  mathematics  beyond  the  level  of  this 
thesis.  We  present  here  a  cursory  summary  of  some  of  the  results  which  are 
consistent  with  the  mathematical  level  of  this  thesis,  and  particularly  those 
which  relate  to  some  of  the  material  we  have  already  discussed.  A  rigorous 
mathematical  derivation  of  the  results  we  present  here  can  be  found  in 
Sparrow  (1982). 

The  Lorenz  equations  were  developed  in  an  attempt  to  model  the  earth's 
atmosphere  to  simulate  weather  patterns  using  a  small  computer.   The 
Lorenz  system  is  defined  as  follows: 
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dx/dt  =  -ox  +  ay 

dy/dt  =  rx  -  y  -xz 

dz/dt  =  xy  -  bz, 
with  a,  r,  and  b  positive  parameters.  This  is  an  example  of  a  continuous 
dynamical  system.  This  system  models  a  flat  fluid  layer  being  heated  from 
below  and  cooled  from  above  (representing  the  Earth's  atmosphere  being 
heated  from  the  ground's  absorption  of  sunlight  and  losing  heat  into  space). 
In  the  resultant  temperature  flow,  x  represents  the  convective  motion,  y 
represents  the  horizontal  temperature  variation,  and  z  represents  the 
vertical  temperature  variation.  The  parameters  a,  r,  and  b  are  related  to  the 
Prandtl  number,  the  Rayleigh  number,  and  the  size  of  the  region  being 
modeled  (see  Figure  4.50). 


COOL  UPPER  BOUNDARY 


zA 
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Figure  4.49  The  model  for  the  Lorenz  equations. 


The  Lorenz  system  is  a  very  crude  model  of  weather  dynamics  and  is  of 
little  practical  value.  Actually  it  has  been  studied  most  extensively  for 
parameter  values  that  are  nowhere  near  those  of  the  Earth's  atmosphere. 
While  the  system  does  have  physical  relevance  to  the  Maxwell-Bloch 
equations  for  lasers,  and  to  convection  problems  in  specially  shaped  regions 
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(usually  toroidal),  they  attract  the  most  attention  because  of  the  wealth  of 
information  they  provide  about  dynamical  systems. 

For  one-dimensional  or  planar  systems  of  differential  equations,  the 
Poincare-Bendixson  theorem  (see  Hirsch,  1974,  p.  248)  guarantees  that  one  can 
completely  classify  the  solution,  as  to  whether  it  approaches  a  fixed  point  or 
limit  cycle,  or  goes  to  infinity  in  a  finite  amount  of  time.  However,  there  is 
no  analogous  theorem  in  three  dimensions  where  many  systems  with 
interesting  behavior  are  being  discovered.  The  Lorenz  system  is  of  great 
mathematical  interest  because  it  possesses  many  of  the  characteristics  of  other 
higher-dimensional  systems.  This  is  not  to  say  it  is  typical,  as  it  has  some 
distinct  characteristics  (for  example,  symmetry)  but  it  does  demonstrate 
characteristics  typical  of  many  general  higher-dimensional  systems. 

Because  the  original  paper  on  this  subject  by  Lorenz  (1963)  fixed  the 
parameter  values  at  a  =  10  and  b  =  8/3  and  investigated  the  system  as  the 
parameter    r  varied,  much  of  the  literature  has  taken  this  same  approach,  as 
we  do  here.  Hence,  we  consider  the  system 

dx/dt  =  10(y-x), 
dy/dt  =  rx  -  y  -  xz, 
dz/dt  =  xy-(8/3)z. 

First  note  the  apparent  simplicity  of  the  system.  There  are  only  two 
nonlinear  terms,  xz  and  xy.  Also,  there  is  a  natural  symmetry  to  the 
equations  given  by  (x,  y,  z,)  — » (-x,  -y,  z).  The  z-axis  is  invariant  because 
points  which  start  on  it  stay  there  and  tend  towards  the  origin.  Moreover, 
when  x  =  0  the  dx/dt  term  carries  the  same  sign  as  y  so  that  all  points 
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which  rotate  about  the  z-axis  do  so  in  a  clockwise  manner  when  viewed  from 
above  z  =  0. 

The  model  has  no  solutions  which  tend  to  infinity;  in  other  words,  there 
is  a  surface  inside  of  which  all  solutions  tend  towards  the  origin  and  herein 
all  solutions  remain.  To  see  this,  consider  the  ellipsoid  f(x,  y,  z)  =  x2/2a  + 
y2/2  +  z2/2  -  (r  +  l)z  -  \i  -  0  for  |i  arbitrarily  large.  We  show  that  the  dot 
product  of  the  velocity  vector  and  the  outward  normal  vector  to  the  ellipsoid 
is  always  negative,  i.e.,  Df  =  (dx/dt)fx  +  (dy/dt)fy  +  (dz/dt)fz  <  0.  Substituting 
into  the  Lorenz  equations  we  obtain  Df  =  a(y  -  x)fx  +  (rx  -  y  -  xz)fy  +  (xy  -  bz)fz. 
Then  substitution  of  the  ellipsoid  partial  derivatives  yields  Df  =  -x2  -  y2  -  bz2  + 
(r  +  l)bz.  For  large  enough  \i  in  the  equation  for  the  ellipsoid,  the  quadratic 
term  in  z  in  the  expression  for  Df  always  dominates  the  linear  term  in  z,  so 
for  this  surface  the  flow  is  always  towards  the  origin.  Hence,  no  trajectory 
originating  a  finite  distance  from  the  origin  will  go  to  infinity. 

We  now  analyze  the  system  for  a  =  10  and  b  =  8/3  as  we  vary  r.  To 
begin/  we  restrict  our  attention  to  small  values  of  r  (i.e.,  r  <  30).  A  quick 
check  shows  the  origin  is  a  fixed  point  for  all  parameter  values,  but  we  would 
like  to  know  whether  it  is  attracting  or  repelling. 

We  introduce  here  the  concepts  of  stable  and  unstable  manifolds.  A 
stable  manifold  of  a  point  p  is  the  set  of  all  points  that  tend  to  p  in  forward 
time  (as  t-><*>).  The  unstable  manifold  of  p  is  the  set  of  points  that  tend  to  p 
in  backward  time  (as  t->-«>).  For  our  purposes,  a  manifold  can  be  thought  of 
as  simply  a  surface  in  phase  space.  These  manifolds  can  be  determined  from 
the  eigenvalues  of  the  linearized  system  near  the  point  p.  The  linearized 
system  has  the  matrix: 
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-a  c  0 
r-z-1  -x 
y   x  -b 

which,  when  evaluated  at  the  origin  yields 


-a  a  0 
r  -1  0 
0  0-b 

The  eigenvalues  of  this  system  are  X\,  X2  =  (l/2)[-a  -  1  ±  ((a  -  l)2  +  4ar)1/2]/  and 
X3  =  -b. 

Since,  for  r  <  1,  we  have  [(a  - 1)2  +  4ot]1/2  <  [(CT .  i)2  +  4a]  1/2  =  a  +  1,  so  all 
three  eigenvalues  of  the  linearized  system  evaluated  at  the  origin  are 
negative.  Hence,  the  origin  is  globally  attracting.  The  phase  portrait  of  this 
condition  is  shown  in  Figure  4.51.  However,  for  r  =  1,  the  eigenvalues 
evaluated  at  the  origin  are  Xi  =  0,  A.2  =  -a  -  1,  and  X3  =  -b.  This  zero 
eigenvalue  is  analogous  to  the  nonhyperbolic  fixed  point  we  encountered  in 
our  study  of  discrete  systems.  We  require  more  theory  to  determine  whether 
the  manifold  associated  with  this  eigenvalue  is  stable  or  unstable.  For  r  >  1, 
the  eigenvalues  are  Xi  >  0,  X2  <  0,  and  X3  <  0.  Since  two  of  these  eigenvalues 
are  negative  and  one  is  positive,  the  origin  has  a  two-dimensional  stable 
manifold  and  a  one-dimensional  unstable  manifold  for  r  >  1. 

As  r  passes  through  1,  not  only  does  the  origin  become  unstable,  but  two 
new  fixed  points  are  introduced  at  (±b(r  -  l)1/2,  ±b(r  -  l)1/2,  r  -  1)  which  we 
denote  C+  and  C".  You  should  recognize  this  as  a  bifurcation.  A  similar  check 
as  we  did  above  of  the  linearized  system  near  these  points  shows  that  they 
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have  complex  eigenvalues.   Furthermore,  for  values  of  r  <  rn  (as  defined 
below),  the  real  parts  are  all  negative.  Hence,  these  points  are  attracting.  The 
phase  portrait  of  the  system  for  1  <  r  <  rn  is  shown  in  Figure  4.52. 


t 


Figure  4.51  The  origin  is  globally  attracting  for  r  <  1. 
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Figure  4.52  Phase  portrait  of  the  Lorenz  system  for  1  <  r  <  ih- 


Numerical  solutions  to  the  Lorenz  equations  indicate  that  for 
1  <  r  <  13.926,  orbits  on  the  unstable  manifold  of  the  origin  tend  directly  to  the 
nearest  attracting  fixed  point  C+  or  C",  as  indicated  in  Figure  4.53.  However, 
for  r  >  13.962,  these  orbits  "cross  over"  and  are  attracted  to  the  other  stable 
point  (see  Figure  4.54). 
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Figure  4.53  Solution  trajectories  for  r  <  13.962. 
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Figure  4.54  Solution  trajectories  for  r  >  13.962. 
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Since  we  know  the  stable  manifold  of  the  origin  is  planar  near  the  origin 
and  includes  the  entire  z-axis,  and  since  trajectories  cannot  cross  each  other, 
the  stable  manifold  must  be  twisted  in  some  strange  way.  What  happens  here 
is  the  stable  and  unstable  manifolds  of  the  origin  merge  and  form  an  orbit 
called  a  homoclinic  orbit.  A  homoclinic  orbit  of  a  point  p  is  a  set  of  points 
that  tend  to  p  in  both  forward  and  backward  time  (see  Figure  4.55).  The 
introduction  of  a  homoclinic  orbit  is  another  example  of  a  bifurcation. 


Figure  4.55  The  homoclinic  orbit  of  the  Lorenz  equations. 

As  r  continues  to  increase,  we  note  that  at  r  =  rH  =  Ma  +  b  +  3)]/(a  -  b  - 1), 
the  real  parts  of  the  complex  eigenvalues  of  the  linearized  system  at  C+  and 
C"  cross  the  imaginary  axis  and  become  positive.  This  is  another  example  of 
a  bifurcation  as  C+  and  C"  become  unstable.  Hence  all  three  fixed  points  are 
now  repelling.  For  a  =  10  and  b  =  8/3,  in  ■  24.74,  numerical  solutions 
indicate  that  for  r  >  i"h,  there  is  an  attractor  (called  the  invariant  set)  to  which 
all  solutions  tend  as  t-»©°. 
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Figure  4.56  shows  this  invariant  set  for  r  =  28,  b  =  8/3,  and  a  =  10  as  it  is 
projected  onto  the  xz-plane.  This  invariant  set  is  a  strange  attractor  and  it 
exhibits  some  interesting  properties.   For  example,  the  trajectory  continues 
forever  within  the  bounds  shown,  yet  never  crosses  itself  or  returns  to  the 
same  point  in  space.   Additionally,  the  dynamics  on  the  attractor  are  believed 
to  be  chaotic,  although  for  continuous  systems  more  theory  is  required  than 
developed  in  this  thesis. 


Figure  4.56  The  Lorenz  attractor.  The  figure  is  from  Holden  (1986,  p.  126). 

There  is  no  closed  form  solution  to  the  Lorenz  equations,  so  most  of  the 
evidence  as  to  their  behavior  has  been  obtained  numerically  and  is 
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conjectured.  However,  the  following  results  for  small  values  of  r  have  been 
verified  computationally  numerous  times,  and  are  widely  accepted  as  the 
system's  true  behavior: 

1.  For  r  <  1,  all  solutions  tend  towards  the  origin. 

2.  For  1  <  r  <  13.926,  all  trajectories  spiral  into  one  of  the  attracting  fixed 
points  C+  and  C". 

3.  For  13.962  <  r  <  24.06  an  invariant  set  appears  in  the  trajectory,  and  some 

solutions  wander  among  the  invariant  set  before  spiraling  into  either  C+ 
or  C".  The  closer  r  gets  to  24.06,  the  longer  some  solutions  stay  near  the 
invariant  set. 

4.  For  r  >  24.06,  some  trajectories  stay  forever  near  the  invariant  set, 
although  for  r  <  24.74,  some  trajectories  eventually  spiral  into  C+  or  C-. 
For  r  >  24.74,  the  fixed  points  C+  and  C"  become  repelling,  and  all 

trajectories  remain  forever  near  the  invariant  set.  These  invariant  sets 
are  similar  to  the  one  shown  in  Figure  4.56,  and  get  closer  to  it  as  r 
increases. 

Using  an  analysis  similar  to  the  Poincare  map,  we  can  see  a  further 
connection  between  continuous  and  discrete  systems.  Considering  the 
homoclinic  orbit  of  the  Lorenz  equations,  we  can  construct  a  small  box  about 
the  origin,  and  analyze  where  orbits  near  the  homoclinic  orbit  penetrate  this 
surface.  This  analysis  shows  variations  on  many  of  the  exotic  structures  we 
studied  for  discrete  systems,  including  horseshoes  and  Cantor  "books"  or 
"fans,"  which  are  families  of  two-dimensional  Cantor  sets  "sewn"  together 
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along  a  one-dimensional  manifold,  or  "spine."    These  intriguing  results  are 
all  presented  in  Sparrow  (1982). 

For  large  values  of  r,  numerical  solutions  to  the  Lorenz  equations  have 
been  obtained  that  exhibit  many  of  the  phenomena  we  studied  earlier  for 
discrete  systems.  For  values  of  r  in  certain  intervals  (called  "windows"), 
stable  periodic  orbits  develop  that  bifurcate  as  r  decreases  through  the 
window.  One  such  window  appears  at  99.542  <  r  <  100.795.  For 
99.98  <  r  <  100.795,  the  orbit  shown  in  Figure  4.57  appears.  In  the  interval 
99.62  <  r  <  99.98,  a  different  orbit  appears  which  has  two  "loops"  that  pass  very 
close  to  each  other  (see  Figure  4.57).  This  is  a  period-doubling  bifurcation  as  r 
decreases  through  the  critical  value  r  =  99.98.  An  entire  sequence  of  period- 
doubling  bifurcations  occur  as  r  decreases  from  100.795  to  99.542.  If  we  let  rn 
be  the  values  at  which  these  bifurcations  occur,  then  evaluation  of  the  ratio 
(rn_i  -rn)/(rn-rn+i)  yields  approximately  4.67  in  the  limit,  which  is  very  close 
to  the  Feigenbaum  constant. 


Figure  4.57  Orbits  for  r  =  100.5  and  r  =  99.65.  The  figure  is  from  Sparrow 
(1982). 
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A  second  window  appears  from  145  <  r  <  166.  The  numerical  solutions 
for  r  =  160  and  r  =  147.5  are  shown  in  Figure  4.58.  Again,  period-doubling  is 
apparent.  A  final  window  occurs  for  214.364  <  r  <  ~.  Period-doubling  is 
again  observed  in  the  solutions  for  r  =  350,  260,  222,  and  216.2,  as  shown  in 
Figure  4.59.  Additionally,  a  symmetric  orbit  is  seen  for  r  =  350. 
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Figure  4.58  Orbits  for  r  =  160  and  r  =  147.5,  showing  period-doubline. 
The  figure  is  from  Sparrow  (1982). 
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Figure  4.59  Orbits  for  r  =  350,  260,  222,  and  216.  The  figures  are  from 
Sparrow  (1982). 
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These  are  just  a  few  examples  of  the  many  observed  phenomena  of  the 
Lorenz  equations.  Additionally,  because  of  the  wide  range  of  parameter 
values,  there  are  even  more  unanswered  questions  about  the  system. 
However,  in  light  of  what  we  studied  for  discrete  dynamical  systems,  these 
results  have  a  direct  analogy  to  the  discrete  phenomena  we  studied  earlier. 

Although  the  Lorenz  equations  reveal  very  little  about  the  weather,  they 
do  give  considerable  information  concerning  continuous  dynamical  systems, 
a  small  amount  of  which  was  discussed  here. 

The  behavior  of  the  system  does  tell  us  that  weather  is  unpredictable  to 
any  degree  of  accuracy  projected  for  any  large  amount  of  time  into  the  future. 
From  what  we  know  about  chaotic  dynamical  systems  (and  weather  is  surely 
chaotic),  even  if  we  were  to  develop  an  accurate  model  and  measure 
atmospheric  conditions  accurately  on  an  arbitrarily  small  grid,  the  sensitive 
dependence  on  initial  conditions  of  the  system  causes  any  computed 
(predicted)  solution  to  stray  arbitrarily  from  the  actual  weather,  given  the 
slightest  reading  error.  Additionally,  small  perturbations  (which  could  never 
be  modeled)  such  as  a  single  person  lighting  a  match,  could  cause  the  whole 
system  to  follow  a  new  orbit.  On  the  other  hand,  if  weather  follows  some 
strange  attractor  on  which  the  dynamics  of  the  system  are  chaotic,  then  not 
only  is  it  unpredictable  (sensitive  dependence  on  initial  conditions),  but  every 
type  of  weather  possible  (topological  transitivity)  could  be  experienced. 
Moreover,  there  will  be  a  dense  period,  providing  some  order  to  the  weather 
allowing  us  to  predict  such  things  as  seasonal  changes. 
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